{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import scipy as sp\n",
    "import pylab as pl\n",
    "\n",
    "import mxnet as mx\n",
    "from mxnet import gluon"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "from gmm_base import *"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "ndims = 2\n",
    "sample_size = int(1e4)\n",
    "num_clusters = 7\n",
    "epochs = 1000"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "gnd_mu_ = np.random.RandomState(0).rand(num_clusters, ndims)\n",
    "gnd_cov_L_ = np.random.RandomState(1).randn(num_clusters, ndims, ndims) * 0.03"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "z = (np.random.RandomState(2).rand(sample_size) * num_clusters).astype(int)\n",
    "x = gnd_mu_[z] + (\n",
    "    gnd_cov_L_[z]\n",
    "    @ np.random.RandomState(3).randn(sample_size, ndims)[:, :, None]\n",
    ").squeeze(axis=-1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_cov(mean, cov, **kw):\n",
    "    vec, v2, _ = np.linalg.svd(cov)\n",
    "    val = v2**0.5\n",
    "    for r in range(len(val)):\n",
    "        handle = pl.plot(\n",
    "            *zip(mean - vec[:, r] * val[r], mean + vec[:, r] * val[r]), **kw\n",
    "        )\n",
    "    return handle"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "def Gaussian_log_pdf(ex, g_mean, g_kR):\n",
    "    model = GMMModel(\n",
    "        ex, num_clusters=1, mu_=g_mean[None, :], kR_=g_kR[None, :, :]\n",
    "    )\n",
    "    return model(mx.nd.array(ex))[0].asnumpy()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "def mixture_full_log_pdf(x, mu_, kR_):\n",
    "    model = GMMModel(x, num_clusters=mu_.shape[0], mu_=mu_, kR_=kR_)\n",
    "    return model(mx.nd.array(x))[0].asnumpy()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "g_mean = x.mean(axis=0)\n",
    "g_cov = (\n",
    "    (x - g_mean[None, :])[:, :, None] @ (x - g_mean[None, :])[:, None, :]\n",
    ").mean(axis=0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(0.2941210180312009,\n",
       " 1.0967109330507165,\n",
       " 0.28322551523446604,\n",
       " 1.0858154302539818)"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAP4AAAD4CAYAAADMz1tMAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nOydd3hU1da4333OzKQTQoBACCQgSC+S0CxgAxvKRVTaFbFhQ6/Xio2fYsXPwrV8KjbQC6gfKCJXvQJSbLQgYGhSA6EkpPdMOfv3x8mcTEsyaYCZeZ+Hh2Sffc7syczae+21VxFSSoIECRJYKKd7AEGCBDn1BAU/SJAAJCj4QYIEIEHBDxIkAAkKfpAgAYjpdL1w69atZVJS0ul6+SBBAoLU1NRsKWUbz/bTJvhJSUls3rz5dL18kCABgRAi3Vd7UNUPEiQACQp+kCABSFDwgwQJQE7bHj/QsdlsZGRkUF5efrqHEqQZEBoaSkJCAmaz2a/+QcE/TWRkZBAVFUVSUhJCiNM9nCB/YaSU5OTkkJGRQefOnf26J6jqnybKy8uJjY0NCn2QBiOEIDY2tk7aY1DwgdT0PN5evY/U9LxT+rpBoQ/SWNT1uxTwqn5qeh6TP1iP1a5hMSksuG0oyYkxp3tYQYI0KQG/4q8/kIPVrqFJsNk11h/IOd1D+ktz7rnn1utaTTz99NO88sor9R1SEB8EvOAP7RKLxaSgCjCbFIZ2iT3dQ/pL8+uvv3q1ORyOaq8FOT0EvOAnJ8aw4LahTBjciX4dopn1zQ4Wbjh8uoflk8a2Rfz73/9m8ODBDBgwgDvuuAOHw0FkZCSPPvooycnJXHrppWzcuJELL7yQLl26sGzZMgDmzZvHmDFjuPzyy+nevTvPPPOM8czIyEgA1qxZw0UXXcSkSZPo27ev2zWAl19+mb59+9K/f39mzJgBwPvvv8+gQYPo378/48aNo7S0tFHeZxBvAn6P7+SLzUewOfQ0ZNsy/gBg0pBOp3NIbjS2LWLXrl18/vnn/PLLL5jNZu6++24WLFhASUkJF154IbNnz2bs2LE8+eSTrFixgp07d3LTTTdxzTXXALBx40bS0tIIDw9n0KBBXHXVVaSkpLi9hrOP5xHTd999x9KlS9mwYQPh4eHk5uYCcO2113L77bcD8OSTT/Lhhx9y77331vs9BqmeoOCj7/OdQu/ku7TjZ5Tg+7JFNETwV61aRWpqKoMGDQKgrKyMtm3bYrFYuPzyywHo27cvISEhmM1m+vbty6FDh4z7R44cSWysvi269tpr+fnnn70Ef/DgwT7PlVeuXMnNN99MeHg4AK1atQIgLS2NJ598kvz8fIqLi7nsssvq/f6C1ExQ8NH3+WZVuAl/7/YtAH2lXX8gh6FdYk+rtd9pi7DZtUaxRUgpuemmm3jxxRfd2l955RXjaEhRFEJCQoyf7Xa70c/z+MjXcVJERES1r+2r/9SpU1m6dCn9+/dn3rx5rFmzpk7vKYj/BPweH/R9/mfThjE4KQYhQADzfjvEwg2HmfzBel79YQ+TP1h/ys/5Pce44LahPDCqe6McOV5yySUsXryYrKwsAHJzc0lP9xnB6ZMVK1aQm5tLWVkZS5cu5bzzzvP73lGjRvHRRx8Ze3inql9UVET79u2x2WwsWLCgDu8mSF0JrviVJCfGMKJ7Wzan56Ghq9OfbzpMhU1D4n7Ud7o0gOTEmEZ7zV69evHcc88xatQoNE3DbDbz9ttv+33/+eefz4033si+ffuYNGmSl5pfE5dffjlbt24lJSUFi8XClVdeyQsvvMCzzz7LkCFDSExMpG/fvhQVFdXnrQXxA3G68uqnpKTIMykRR2p6Hku2ZLA4NQOHQ0NVFTRNw67p182q4Jlr+vD0sjRsDolZFSyaNqzegrhr1y569uxZp3tKKuyUVNiJCDEREXL65ux58+axefNm3nrrrdM2hiDe+PpOCSFSpZRes3JwxcfdYm5SBBMGdyKrqIKVOzONPpomWb0nC2ulHcDqkCzZknHKVv2SCjsHs0uM/XHn1hFewn+mTAxBznxq3eMLIT4SQmQJIdKquS6EEG8IIfYJIbYLIQY2/jCbFleLuUOTZBVV8OPuLFx1ISnhYHaJ233ZRRWnbIwlFXaklEh041hJhd3r+sHsEjILyzmYXeJ13dknq7Dc572+2qtj6tSpwdX+L44/xr15wOU1XL8C6Fb5bxrwTsOHdWqJCbegCIECqIrgx91ZOLQqsReAxazQubW7lbp1VMgpG2NEiAkhBALdgu65oueVWtF8TAxOoc4prvA5MfgzYQRpftSqD0op1wkhkmroMgb4ROrGgvVCiJZCiPZSyuONNMYmJTU9j1nLd+DQJELAgI4t2XSoynqvCJgwuBPjBiYAsHZPlrHHd7bV5zUd5TZDLXelJnW9ZbgZAbQMtxjXSirs5JdaySu1VXWsnBhctwcgkE4dpnJicPZxahKu7UGaN43xCXcAjrj8nlHZ9pcQ/PUHcgzLvZSw8ZD7kZ0Axg1MMPbyi6YNa5BV32lPeOuKOA5ml7jt1Usq7Bxw2cd3qbzmub8PNauUVNhRFcHxgnI0DwNtmFkxnucUaoFE1xekm8bg1CSonBxUJRgqHAg0huD7+qb4PCoQQkxD3w7QqdOZ4RUXE27xPdhKHBJun7+JlKRW3DHirAYfqTntCa4quVMI80utOE9ZpJTkl1p9rsrH8stx/ol9jb3U6uBAdgnx0aHgFGohiI8OpdzmcLsnIsRE++hQjuWXI5EcLygn1KwGV/1mTmM48GQAHV1+TwCO+eoopZwrpUyRUqa0aeOV4/+0kFdqrbVPbqmNH3Zmcv27vzY4gMfpgedrr+4pxM7f9VVYGDOsrFTaPfu7zsBSSrKLvd9bXqmN3BIrB06WkFOsGyd1e0bVhOPvPv/CCy/0WRuhunZ/cA3BnTlzJitXrgT0OgzZ2dl1ft7DDz9Mjx496NevH2PHjiU/P99nv6SkJPr27cuAAQPcfBKeeuop+vXrx4ABAxg1ahTHjulf7d27dzNs2DBCQkL+kiHDjSH4y4Apldb9oUDBX2F/74x021CH+HtNwlNfpzXIg8/pgdcizOR1JBdmVg3hdbq0Hi8oM1ZjqEaVqsTzWoXd4aZBZBVVGNsCia45ODWOmgyHp4tZs2Zx6aWXNugZI0eOJC0tje3bt3P22Wd7uSi7snr1arZu3eo2aT388MNs376drVu3Mnr0aGbNmgXo8QVvvPEGDz30UIPGd7qo9RMWQiwCLgRaCyEygP8HmAGklO8C3wJXAvuAUuDmphpsY+F6bq/VJEk+cGiS99buZ+4U/z3VPElOjGFX6Ymq/XUTsu1I1SRlc2hu196b8zL/XbaYpMROtIyJpXe//jz88MNcddmlDBkyhNWrV5Ofn8+HH37IBRdcQFlZGTfffDM7d+6kZ8+elJWV1fr6ixYt4oUXXkBKyVVXXcXs2bMB+PDDD5k9ezbx8fF069aNkJAQryPCqVOnMnr0aK677jqjraysjLFjxzJu3Dgjkq8mRo0aZfw8dOhQFi9eXOs9rrRo0cL4uaSkxPi82rZtS9u2bfnPf/5Tp+edKfhj1Z9Yy3UJ3NNoIzoFuJ7b14dVuzJJTc/7y6boEkDatt9Z+d0yVv+yAZOQnDtkEOcMHGis9na7nY0bN/Ltt9/yzDPPsHLlSt555x3Cw8PZvn0727dvZ+DAml02jh07xqOPPkpqaioxMTGMGjWKpUuXMnjwYJ599lm2bNlCVFQUF198Mf3796913MXFxUyYMIEpU6YwZcoUAC644AKfrr2vvPKKl7bw0UcfMX78eN9/EyEYNWoUQgjuuOMOpk2bZlx74okn+OSTT4iOjmb16tW1jvOvwJmh051inPvscptWe2cfOCQNDos9XQghaB1p4fdN67lw1JXkWXX7wfkXX0ZuidXY31977bUAJCcnG+G469at47777gOgX79+9OvXr8bX2rRpExdeeCFOe87kyZNZt24dACNGjDDCca+//nr+/PPPWsc+ZswYHnnkESZPnmy0/fTTT3697+effx6TyeR2ryu//PIL8fHxZGVlMXLkSHr06MHw4cONe59//nlefPFF3nrrLbfEI39VAiI6z1fmmn4douv9PCFolBRdp8NZRkqJpsnK4zv9d18egc5wXFVVawzHre216tJeG+eddx7fffed2/0XXHABAwYM8PrnNAoCzJ8/n+XLl7NgwYJqxx8fHw/oKvzYsWPZuHGjV59JkyaxZMmSeo39TKPZr/iemWtmju7N09/swGqv32oPIBoprqmkws72I3k1GuyagqJyO+cMHsasGfdz2/QHcDjs/PTjD1w36aYaDXvDhw9nwYIFXHTRRYbBrCaGDBnCP/7xD7Kzs4mJiWHRokXce++9pKSk8M9//pO8vDyioqJYsmSJkZ6rJmbNmsWzzz7L3XffzTvv6A6ita3433//PbNnz2bt2rVG4g9PSkpK0DSNqKgoSkpK+OGHH5g5cyYAe/fupVu3bgAsW7aMHj161DrOvwLNXvA9M9d8l3YcWwOEHkCDRgnQOV3OMlaHRp/+A7niytFMumI4HTt1Ijk5maT2rWsU/Lvuuoubb77ZON4aPHhwja/Tvn17XnzxRS666CKklFx55ZWMGTMGgMcff5whQ4YQHx9Pr169iI72TwObM2cOt9xyC4888ggvv/xyrf2nT59ORUUFI0eOBHQD37vvvsuxY8e47bbb+Pbbb8nMzGTs2LGAbtuYNGmSkYVoxowZ7NmzB0VRSExM5N133wXgxIkTpKSkUFhYiKIozJkzh507d7oZA89kmn1Ybmp6HhPfX4/NrmFSBRd2b8uaSrfbhiAEjOwZZzj11HVMZVnpRLdPOuWrvSuqvYLeSXGUlpYyfPhw5s6dW6vBrrEoLi4mMjISu93O2LFjueWWWwzhC1I/6hKWGxB7fOd+1uaQrNiZiaZJn+6GdX3kDzszuaGOTj2p6XmMf+9XyirdhE8nMx6YzoABAxg4cCDjxo07ZUIPuqPOgAED6NOnD507d+Zvf/vbKXvtIAGi6ts16SZkDVzs3XBImPl1Gt3bRdW68qem5zHrmx00cKfhhlZagGYtx9Qyrs73zv34E9pHhzXeYOrAX9HbrTnR7Fd859FdU6Jp0qsCj+dJgtPIuC2joMZnSVsFjqLavQml5sCWdwx74Um08iK0ssI6jzun2BoMww1Qmv2K73SRnfXNjlqFrr4oiiAm3GJk5I0JtzBr+Q4qbBqqIpg1pg87jhVQUYvfgKM4F0exnnhSmEJQwiJ99tPKS7AXZoHmMNrshdmYLWEI1b/66OAdJBQkcAiITzw5MYaZV/dm/Hu/1qhmn/x6NuWHfgfVhFDNxF5+H2Gdz6n1+Q5N8vQ3O0BK7JX2A+d2wq5JnvzqDxA1+9nbC7LcVm17YRZmSyhCdfmIpIa9MLva1V3abbUKvqtt40zyyw9yagmYTz05MYaLe8Txg0sePU+0ihK08mKXBke1fV2RUKNfgObsVANqREu0sqKqjlLDXpCJuVUH/VdrObaCTHDYvO4V5lBMLeP8Wu1bR4XoSUeoSugRzNUXeDT7Pb6T1PQ8Vu/JqrGP9BQqk/9qc0MRJgtqlLs3oLSW4SjJw1GUgy03w6fQq5GxmGMT/Fbxs4sqyC2xklti4/DRY9wwfiLdz+7KhecP5dxzz2XRF3ULYqkPmzdvNlx/G8qZFAKcm5vLyJEj6datGyNHjiQvz3cU5yOPPELv3r3p2bMn9913n+GJaLVamTZtGmeffTY9evQwvAQrKioYP348Xbt2ZciQIW4VjepLwAi+rzJZXjjcDV1uavYpQI1oibC4W9kdRTk4Sry/QMJkwRzbETWybj4Ezr+AJjUmXH8d/QcN5dtftvLZt2uY/fYH7Dvof1GN+pKSksIbb7zR5K9TFxojBPill17ikksuYe/evVxyySW89NJLXn1+/fVXfvnlF7Zv305aWhqbNm1i7dq1gB4T0LZtW/7880927tzJiBEjAD2SMSYmhn379vHPf/6TRx99tEHjhABR9VPT81hTy2oP3it+XQxl9aXfB4lN+vztt/kW5I2/rMNsNjN+yi2UWvUtTXxCJ26ZdhcAhw4d4sYbb6SkRM8s/NZbb3HuueeyZs0aXnnlFZYvXw7onnEpKSlMnTqVGTNmsGzZMkwmE6NGjeKVV17h//7v/3jmmWdQVZXo6GjWrVvn9oyNGzdy//33U1ZWRlhYGB9//DHdu3dn3rx5LFu2jNLSUvbv38/YsWNr9dQ73SHAX3/9tVH266abbjKKj7oihKC8vByrVc+2ZLPZiIvTj2I/+ugjdu/eDegly1q3bm089+mnnwbguuuuY/r06dWWIfOXZi/4zmM0fyLxZANW/K5tIjiQXVLvUN9TiUCw/8/d9Ozbn5hwC2W2cuOL1DLcAujBKitWrCA0NJS9e/cyceLEGlXq3NxcvvrqK3bv3o0Qwsh0M2vWLP773//SoUMHn9lvevTowbp16zCZTKxcuZLHH3/cUHG3bt3K77//TkhICN27d+fee++lY8eOXs+AMyMEODMzk/bt2wO6u7KzPJkrw4YN46KLLqJ9+/ZIKZk+fTo9e/Y0/jZPPfUUa9as4ayzzuKtt94iLi6Oo0ePGu/bZDIRHR1NTk6OMTHUh2Yv+E5ffX+oz4ov0KvsHM4r+2sIfWXuvahQE9FhZmIjQwg1q9w7fTqbNvxKaEgImzZtwmazMX36dLZu3YqqqrWGzbZo0YLQ0FBuu+02rrrqKkaPHg3oEXVTp07lhhtuMEJ9XSkoKOCmm25i7969CCGw2ao+g0suucTw4e/Vqxfp6enVCv6ZFAJcE/v27WPXrl1kZGQAeoagdevW0atXLzIyMjjvvPN47bXXeO2113jooYf49NNPfUYzNjSBS7Pf4zsdePx5o54rPrUI/t8GxPPQZd25PqUjdkcjuuM1JVJSbnPQvUcvtv7+O6An3Pzo/XdZ/eOPnDx5EoDXX3+duLg4tm3bxubNm7Fa9fx9JpMJTat6r+Xl5Ub7xo0bGTduHEuXLjWCXN59912ee+45jhw5woABA8jJcXdOeuqpp4xov2+++cZ4HlSFBoN3eLD32zr9IcBxcXEcP65nnTt+/Dht27b1eu5XX33F0KFDiYyMJDIykiuuuIL169cTGxtLeHi4Ea9w/fXXs2XLFgASEhI4ckRPZG232ykoKDAmsvrS7Ff85MQYZo7uzcyv0/Q49JrwWvGr//NEWlTmTNDP+BduOIxSmc3WZFJwaBr+zgPbpuyuPKbz/lILkwVpd0+YqYRFYYquu3uucMm2m1tq4+xzhlJYXMqcN9/i/nunAxjVa0FfiRMSElAUhfnz5+Nw6HaAxMREdu7cSUVFBeXl5axatYrzzz+f4uJiSktLufLKKxk6dChdu3YFYP/+/QwZMoQhQ4bwzTffGF9g19fp0EE/spw3b16d35eTMyEE+JprrmH+/PnMmDGD+fPnG5GIrnTq1In333+fxx57DCkla9eu5f7770cIwdVXX82aNWu4+OKLWbVqFb169XJ77rBhw1i8eDEXX3xxcMX3B2eVmdrw3uNXv+KX2hykpucZBTmc8QCJMWF+Cr3EUZSNLfeoD6EXqFGtMbfuhBrR0u2KVlbk7mvgBwKIjw4lrkUoMeFmYwKY88G/WbNmLZ07d2bw4MHcdNNNhjHq7rvvZv78+QwdOpQ///zTqHXfsWNHbrjhBvr168fkyZM55xx98isqKmL06NH069ePESNG8PrrrwN6ssq+ffvSp08fhg8f7rW/fuSRR3jsscc477zzjMmlPriGAPfv35+BAwcyZswYOnToYIQAX3rppXUOAS4vL+eRRx7xq/+MGTNYsWIF3bp1Y8WKFcyYMQPQjy9vu+02QDfOnXXWWfTt25f+/fvTv39/rr76agBmz57N008/Tb9+/fj000959dVXAbj11lvJycmha9euvPbaaz5PC+qKX2G5QojLgX8BKvCBlPIlj+uJwEdAGyAX+LuUMqOmZ57KarlOA5+tMp99dQv/4devR1qrEkh2vP8LlBDfyRsALujWmo6twuucctt68hAfje9O61Ytva4JkwVTy3YIk25kQ0psOUfcV36hYG7dqU7Gx1YRFhJiwv0qvtncCJQQ4EYNyxVCqMDb6DXyegEThRC9PLq9gl5Gqx8wC6g+h/FpwOmvP2FwJ2rKfeFt3KtZIH7am81nG+sm9FJKspe/Cpq3aq9EtMQc27FK6AGEwNSyHa7OtsJcWSijDjh7R4Toab3jWoQGhNBDMATYF/586oOBfVLKAwBCiM/Q6+XtdOnTC/hn5c+rgaWNOcjGIDkxhvUHcozVXjrsWDP3I1QzQjUjFdVL5ZaKUmvcfl0t+UIIYi/38FpTTJhaxqFYfIfIOr36HMU5qFGtUcPrli/Q9ZgOCDjX3GAIsDf+fPq+auMN8eizDRiHvh0YC0QJIWKllG4m3NNdQmtol1gURaA5JFp5MSc+fbDG/kf+52+gqIjKoJ0Od81DsYQ2eBwh7buBRd9CKKFRmFq0AaVm5UuNaIkSGlEvp6KYcHNACXogUtfTC3+Me/7UxnsIGCGE+B0YARwFvHTZ011Ca8+JIuyVbrtefvnVoTmQtgq08uJGdeE9UmZGhrZAjW5bq9A7qY/QK0IPGQ7SfJFSkpOTQ2io/4uSP9/kWmvjSSmPAdcCCCEigXFSyqYJfm8An2+q2o97ndnXigBFbbSxvLkxn3uFILFlkUtVvPpjUgQOl0xDAggxK7QINXO4MCAObwKa0NBQEhL8L9vuj+BvAroJITqjr+QTgEmuHYQQrYFcKaUGPIZu4T8jcE2OkXasai4SioKlXTekw4Z02JH2ChyFJ6t/kNq45a4KKzSeX+d/3b66YlIEn98xjP5/waIfQZoef0po2YUQ04H/oh/nfSSl3CGEmAVsllIuQ6+t96IQQgLrOENKarnm1FeEcDtfN0XH0f6m143f7QWZHH33VuN3tUUbOtz5ka7qO+xIH1b4MxVTZdafv2KlnyCnBr82rVLKb9GLY7q2zXT5eTHQ9IHcdcS9Rp6efKI6E4i0e/vpCyEqs/GcWYYxVQGknuVHFSAUgdQkqiK4PqUj1w5MCAp9kBo5s77RjYzTT99m11BVBbtDozrjp+eKfqYJuyttIkPILKwwfr8hpSMdWoYxtEtsUOCD+MWZ++1uBJyOO+sP5HA0v4xFNXjY+Vrx/UER7mf5vdpHsfO4d+hmY3LCRejNJoVxwRU+SB1p1oIPuvAnJ8aQmp7H4s1HsFaThcc7Mq/2P40CPPc3PeDju7TjXNGnPd3bRfHkV3+w60TjCL8ioG+HaLZnFPjcpkwdlhQU+iB1JqDOea5P6YhFrcYyX8/sO3mlViYN6cSntw6he7soJr6/vtGEXqBPLDOv7k2I2Xdo8fc7TrhVAQ4SxB+a/YoP7tb96lxs65WEw6Ncdl2SftTG4KQYHr2ip/HcmaN7k1dqZW9mEUu3VrlRHMopZfIH61lw29Dgyh/EbwJixXe37vumXmm3PJQHTw85BbhzeJc6V8Ud1SuOL+48F4DJH6zn1R/2MGv5DoZ2iWXOhHN4YWxfkmKrogZtds2rkk+QIDUREIJfVGarNZimPiu+psGclX8aqvaOY+7Oipf2imPGlT0ZP6ijMUcoQP+EaCwmBVXox3GetInSM8+4TlhWm2a81qQhnXj1hgGEmvVnmE2Km+YRJEhtNHtVPzU9jw9+PlhrPy/ffdXMMGUHB7T2ZOI7zZEEft6bzaZDuSy4bahP49vbq/fRJz6aELN+rGg2Kcy8ujeA4VE48+s/jAo/ZlVw7UDd9dJ5HGm1aWjAL/uqXsv1xCJ4jBekrjR7wXdWy60VD1X/qoQi5plnc0S2Zbz1KXLwHQorqVK1xw1MYPHmI9gcElUVrNmTxcpdmVhMirFHdxVS5//d20WxZEsGAtycb5zCPWfln/yyLxtNVr2W87QiKPBB6kOzF/yhXWIxKaJW4Xdd8a/rZWLB0P2YhaSrOMa/LS8w0fok+UQZfRShR75JKQ1VOzkxhkXThhl+A59tPGwIq1PonXtxV4GtSYCTE2O4/9Kz2XQo19AYgmp9kIbSbAXfGZwztEsss8b04cmv/qAme7sa2YqQTv3AYWNicjZmpcpJpqdyhE8sL/F36+MUoueeq2kVd/oNfLklA6tNQwhBUZnNOFmwmJQ6WeGDan2QxsavnHtNQVPm3EtNz2Pi3N+wOSRmVbBo2jAA3l27n9/2Z1NcUXNSRws25ppf40J1m/tztW7caH2MUkK5c3gXZlzZs5on6CzccFjP7islihBoUqJV+tc/MKo791zUtWFvNEiQWqh3zr2/Iku2ZGB16LHpVodkyZYMkhNjeH9KCo9f6Zku0BsrZu6w/ZNfHe59k5W9fGh+hVAq+O1ADtM+2cyYt36uNtmmM7uvJkHTdOEPWuGDnAk0S1Xf84Qsu0hX21PT88grtdI2ykJWkdX7RhcqsHCb7SE+ES+RolRVXhmm7mQurzHt6ENsy9CP/LZl/AHApCHu6cRcg4TM1WwNggQ5HTRbVX9CpaoPYFEFt5zXmQ9+Puifhd+FKEr5t+UF+isH3NpXOJK5y/YP7JVz5wXdWvPprZ6pCN1tDUFhD3KqCShVPzkxhutTqpxm7A7J3J8O1FnoAYoIZ4p1Brs099V8pJrKHPPbqOj2giv6tK92LPdc1NUw+L29ep9RiMP5c5Agp5pmqeoDjBuYwJdbMrDZdau6owEVLQuI5O/Wx/jM8hzdlKNG+2h1AxHhERy/6HUvNd8T13gBkyJACOyOulv4gwRpDJrlig9VR2APjOrOrDF9CDHXniO/JnKIZrL1cQ5p7nXrLqr4kUknX6faDB+VLNmSQYVNd7+1OSS2Sldcf/3sPTWE2jSGoEYRpCb8WvH9KKHVCZgPtKzsM6MyXddpxdUxpnu7KNYfyKGozMb8Xw9RVo8oulJLayZZn+CLkFkkiOyqC6nzwBQKl79kVLhx3dsDLE7NMFx6VQUURcHhqN4hx/N+Vx+AmaN7M2v5jmp9Aly1i6BGEcQXtQq+SwmtkeiptjcJIZZJKV0r6TwJfCGlfHzKEgMAACAASURBVKeyvNa3QFITjLfeOL/4E99fX+/QWVUVHKdS+C2zaCdcVtMN7+rCf+nTpB7OZ/IH66mwaaiK4OIebY0y2gIYP6gT1w5MqNbo5ym44wYmGME6NrvGd2nHq4J37Hrwzv2Xnm08xzW4x9XFN0gQJ/6o+kYJLSmlFXCW0HJFAi0qf47GI+/+mcKSLRkNipcvKLMjgcMyjsnWx8mWHv77v8yBtbNZfyCHCpteoNOuSVbtysSk6pF0IWbFCMKpDuf9zqg8CUY0n9mkcEWf9lhMipH266e92Ux8f72h1juPEXVTQrCgRhBvGquE1tPAD0KIe4EI4FJfDzqdJbRS0/P4v81Hau/oJ/tlB260PcZX4S8SancJx13zImNSJK8r/Y1TBAlceHYbymwOw/o/+YP1lNs0FAHTLujCyN7tDA0gJtxibAs0oE98NOMqNQRnfYDh3dpwILuEfVl6yWyrXeO9tfvp37ElQ7vEMnN0b2Z+nYZDk8xavoPu7aJqXfWDR4+Bgz+C708JrYnAPCnlq0KIYcCnQog+lQU2qm6Sci4wF/Rz/PoMuL6sP5BjlM9qLPbSiVfiXuK+jAdpIUqN9oTNL/FZ/xlM2NYfh6aXo/5xTxaaJtl0KJdrByZQbtP/NJqEd9cdYO5PB9Ck/seODKmq2COAtGMFxqnBxLm/GXkDPfN7rNqVaUQDjhuYgCalW/RgTcJ8JtkFghNQ0+OPql9rCS3gVuALACnlb0Ao0LoxBthYDO0Si6m6fHt+0iaySmXWc/QLPtwfzVTrI5TIELe+KbteYsHA3aiVpa3sDmnsyXccLfCqcm1U8QWKXGIJJPD5piMs3HCY9QdyDKck5z2iciyq0Ps69/We24PaXIR92QVOB84J6NUf9jD5g/XBU4kmwh/BN0poCSEs6CW0lnn0OQxcAiCE6Iku+DXUozr1eDr11JWW4Wb+ObK7kfVGcalVt0Weza3WRyiX7ll7BqfN4mp+cmvTJGzPKKjTOapDk8z8Oo2YcAtmj8lLAqoiuP2CLm77+j7x0cZxpj+rt9Mu4DpRnI4jwTNlAmru1Pr9k1LaAWcJrV3o1vsdQohZQohrKrs9CNwuhNgGLAKmytPlC1wD4wYmYDbVz3Uhv9TGU1+nMXVYEhMGd8J1t2NSQEs8lztsD1Ihq3ZPAskr5ne4Ulnv9iyJfuw/KCmm+qy/HmiaJK/UyqJpwxjZK85tAnNokqgwM1OHJQG6QXHW8h0ARg6AhRsO+xRip3ADbhMFcFpWXl8TUJDGp1n66tfE41/9waINh92qytblL6AKuKRnHD/szDTaBifFsOVwPnZNcomSyrvmOZhFlbpukyp32e5npZbs9iyzKtxU95owqYKk2AhahZs5lFPiFWR05/AuXrEIXdtEcDivzFD9nRV0nRpATfv6t1fv49Uf9pyWMOLgHr/xCChf/ZroEx+NqggUAaFmhTED4ut0v0PCxoPu6mduidUQuFVaMv+w3YNDVq3JZuHgbfO/GK5UxfcL8FvoW0dasDsk+7KK2Xgoz2dk4dKtR73ckvedLMFaKfSgT3BWF/W5JrX6dK68rvENQZqGgBL81PQ8Zi3fgaMyNn7qsCS+33Gizs/JL3PPz9cqwv2c/FttKA/a7kJzEf4QYWeu+TWGKrrfk79ahgCyi2sOIQa9rJY/z1SEMIS4JuF2dXkOev41P5ptkI4vnCucvseW7Dhe2CgFMDYf8t7/LtXOJ8RuY7b5faMtVNj40Pw/3Gh9jC3ybL+eXdeNWP+EaP7IKHBLM2ak9lYEF/Voa7S7pvSKCbd45QMMJvNsvgSU4HtWzw0zq3oiToesMR9fbVR3739MlxJqs/KMeb7RFiEqmGeZzSTrE6TJLg14Vd9Eh5m9DBdCwKU941izJ4tVuzJZuyfLq5z2mXKGH+TUEFCqfnJiDDNH96ZPh2g0KVm5KxOEYMKQTrwwti8h9bT4V0dxhYP5jst40TbRrb2FKONTy0v0EFUpuxrLyLpub7ZX8RBN6sU+rE5fAodk4YbDhrU+eIQWeATUiu/c4zv96AEcDo0OLcM4nFNCRSPVvfPkPcfVhGLln+YlRluMKObflhcYb32K/bID+es+oeSPlajRbTFFxxHZbxRhSQMabQxH88vdfnc19HmmCAseoTV/AmrFd93jO3EGsSzderTa+xqDfzmu5R371W5trUUhCywv0ElkYi/IxFGSh/XYHkp3rcNR2PT+T5qEbUfyAYKGvAAjoATf1YptUnSPN01WOrs0uTuDYLZ9Ah/bL3NrbSfyWGh5nvZ2dy9oU7R7wo+m4oedmUz+QHcw8pUiDCD14MlgUo9mRkCp+q5W7GP5ZSyqrHRTYdMoEvbaH4DuzHJe19as25tde2cvBM/YpxCCjUmmH43WBJHN91fA8OOC48X6DKRGt63uIY2O1SWIx9Op5/MeP1P650+8aX2AN00hQY2gmRBQKz5UOYf0jo82jrkkUGKtuciGE4eEzQ1a+QRP2G9hieN8t9auMbBqSjhtwgUIBVOLNg14jbrjjNmvMvRJ7pOL6L/3LS4Q23jT9AbSbg0a/poJAbXiA0Zpq//bfIT6RumW+jlJVIdE4RHbHYRgY7S6wWjv2UZlxY3hXPplOEJRa3hC46JJjJh9fTskeER+yi2m74w+I9VUXhQfkthlfg1P+usTKO7CASX4TjXW1ap/unCgcr/tHkKwM1JNNdr7t1P5/gbJVEopIrzJXl8AibHhpOeUusXs3zOiCz/1Wk6bPd+59S9XIuh15X30aMbCcCblJGhqAkrV97TqC/RiGxEWFUdJHra8Y0it/qu5Lf8EjrJCv8/k7Zi4x3Yfq/Ldc/Int7byseVlwimv5s6GI4HLe7cjxKygoJ9utApTYdl02uxZ4NY3X0bwd9vjlMTpQUbNNYNvIPkzBNSK7+m5d11yAn3io3l9xR6O/bGK/LXzQDFhjmlPVMoYogZcXqfnn/zqBWxZBxCWMEzRcZgqz+RNLfT/nWf0SmgkojIThxUzkzf2YlG3I1zUuerjSFH+5EPzK0y1PUIFDc+ZF25RvbYoUWFmI0WXIm20+PYuUH5z65MtW/B36+PspZMhCM11VQwkf4aAEnzQY/Jl5f9QlfvOlpOhd9Ds2HKOIO21B8Z44ijQQ3WltQzbyUPYTh7y2c9zYig/ksbVm0v54cZwzu1Y9ZE46/TdbnsQK2afz/KXCItKmdVhaDuWyi/2+gM5mKSVf5ne5DLFPUw6U7ZkkvUJDtLBrb9zVaywaTz4xVamDT+r1oIijUVT7sEDqRx5wAi+a8lqZ04655cYwJ6b4dbf3KpDnZ6vlRejVZT41dfXxGADrlhQyqopEaTEVxn2RqjbeT3vMW4/Ngai41Gj41BCIgyNwV9OukT4DU6KoWW4hSVbMugfF8J7ltcZIba69T8qWzPJ+jiHZTv6JUQz8+rehiBYTIphJzmUU8rjX/kuGtrYnIo9eKAEJgWE4Kem5zHz6zQjZr7CpvHllgx6x0ejSd1P3uYp+LEdfT2qWhzlxZhjO2EvyETaK+o1zsIKuOzfpay+KZx+cVXCf1XrY7x+4g0mzy/DIUFYwqu2Eca/toQk9EINj67hFXQ2HcpDAmGUM9ryKucqO9yup8u2TKp4gqPoR4q7ThQZ15yr4oNfbOVQTlWC0e/Sjje54AfrBTQeASH46w/kuCWpkMDnm4/Qu72eFlsrK0QrLzauC1MIaou65Qo1t2xH/G3/i5QSrbQAe0Em9oIs7IWV/xdk4ijIwl6QVePEkFsmufSTUtZODadnmyrhH9/HTIVDMnVpOdJa6nMr0fb6ZwjrkkxtSCCSUj62vMwglxLgAPu0eCZbHyeTVkabw+EuZMmJMUwbfpax0kP1RUMbk0Dagzc1jVVC63Xgospfw4G2UsqWjTnQhjC0SywhZsVIaQ16Bd1tGbrge672plbxCFG/Aw8hBGpES9SIloTEd/e67jkxVBzfQ9GmpW59TpZKLvmklHU3R9C1VdU4pvS3UGaHO5f7tvb76+bbkmLmWV5igEfp7+OhZ3FTycOcFC3cXJhV1VvInKv7d2nHuaJP+1Oyxw+kPXhT0ygltKSU/3Tpfy9wThOMtd44vzBLtmSwODXDyEHnxDDsVWJuVXOlm4bgOTEoljA3wQ/p2Jc2Yx7BXpDF9cf2sSz8SzqGlhnX70i2UG6X3P+9t9ag+uHt14pCFkfMpovjoFv7dq0zy856g3s6deK7tOP8vDfbyNN3XXKCTyGbNKSTm8CfCueXQNmDNzX+rPhGCS0AIYSzhNbOavpPBP5f4wyv8XB+YZzlsz/fdBhnFK491z0yr66GvYZgL8xy+90UHYcaEYMaEUMe3fm7SOYLOYs4kW/0+ceQEJSzL2bWoT7GNkIrL0Exh3g+3o025LHQ8gJdHO7vN1Xrxi32GRRtKcCyfQczR/dm06FcQ6UeV0vJLzg1hrdA8ao7FTRWCS0AhBCJQGfgx2qun7YSWk6cE8DJogojU663Ya/pVnxP7AWZbr+bPIJz0mU7Jlmf4HPLs7QWhUb7vTG/YItszxuO6/16nXiyWWB5ns6K++ttVfsyI2QGBQX6KYHNrpFXaq2zSu1peFuyJaNRhTSQvOpOBf5sZP0poeVkArBYSunT/U1KOVdKmSKlTGnT5tQGobiSmp7H6j1VK62nqm9qQlXfE3uB94rvyX7Zgb9bHydPRrq1P2BezDT1m1pfo6PI5IuQWV5Cv87Rl0mlD7C3oOoj9rWf9wfXkGdVESxOzfDKyV9fj7/U9DzmrPwzYLzqTgX+rPj+lNByMgG4p6GDampc6+hJhw17vnumXU9Vv6659+tCbSu+k92yEzdaZ7DQ8jwtRNWe/3HzIiqwMN9xmc/7uohjLLQ8717SG1jhSOYe231ejkHXJeuT3sT31xuq/qLba19dqwt5dhXS+qzYnvEVSuXEciy/jNT0vOCqX08aq4QWQojuQAzwm+e1M42hXWKNijr2vBPgUttTjYxFsYS59W/KgB5/VnwnabILU62PetXpe8Y8n/Hqaq/+3cVhPrfM8hL65Y4h3GX7h5fQK0L3aHxv7X4jpsFq130e/MEZ8nztwASvtN319YN3ja9QgL4dokEIFm08HKyt1wAaq4QW6Ea9z87E0lmuOA1Et5ybBPja3586w55mK0crrTLaIRTUqJrV7C3ybG61PexVp+9F0wf8TfnZ+L23OMhnludo42IXAPjKcQH/sE3H7kPZu6SnPums2u0+Gbl+oNWp667tvnLy17dAh+t9FrNCnw7R2B1Blb+h+HWOL6X8FvjWo22mx+9PN96wmgZXA5ETm4dF39SqalfTlCo+gKPAPa+e2qKNX3H467VeTLM9wPvmVwmpzBykCMmr5neosJk5IVsx3zLbrXQ3wCL7RSzt8BCO9AKvZ6oC7hxxFusP5LhFF6qiKq6hOgNbde2uabyGdomts8HQOUnPHN2bvFKrMVks2ZIRdOJpIAHhuefEVd104n2GX7Xi92gXxe4TRW7Crwg4v96pt9zxd3/vi3Vaf+6x/YN3XOr0qULypvlN7KiECptb/4/tl/GMfQrq4QKjdLcTIeDZv/X18sVH4CagS7ZkGHttV5fZ6lxpXffnqiKYNaaP3/X3arLiu04ggDGxBPf7/hNQ8fiuaqPTju0VnONylLc7s8hrxb+mf3yN6qVSh9iZ0MR+tL/lLdqMm0nMpXcQ1b9uYcArfdTpMwnNS+jfsV/NM/YpgMAhoU98C/cxA93bRQFVtQcURSAlbDyUx8S5v7Fww2EWp2YYfw9VEUYp7WP5ZZgU4aXGrz+QY0wU9spS3/7uyWuyCThtCXB6Kvo2BwJqxXe1PO/NLOKr34967/Erj/IEeilrTw5ml9RY7NKzmEVNCJMFS5skLG2S/L/Jg2+1oYTYbLxqfsfnpPO6bRz/clyL66lsdJgZRVSNVUrcCmkeyy9Dc3kjNofku7Tj2B36FkkA16d0ZM+JImZ+nYZDk5hNCuMHd2ScS3WeoV1iURVhBEdpUvodWOOPX34waKf+BJTgQ5UDzxNf/VFjcE518hvXIhSzWoi1vgn7moByLEgfFgmrVPlJ64unK4brNkWgG81iwi2Gam1SBCaXEt5mVdC7fQt+258D6GHNUSEmnlz6hzF52O0aJ4sqmLPyT3q3b0FUmJmhXWKZNaaPWzi0U0uoba/vj19+MGin/ojTZYRPSUmRmzdvrr1jE3H7J5v55ocfyVzwqNFmbtuZ+JvfBEBV9FM+IQAh0CpXtUW3DwXgpe926dl2ZcMMgJEhKsUV9U/3NUb5mdcs76JWU8GvUIYz0foEO2Rnr2sCOKtNBEMqBcZ59i6Akb3iaBMVwskivQrvj7szcejbfgYlxZB6ON+rLLcnoWZ9bw647ckb0wMv6MZbM0KIVCllimd7wK34oH9ZVuzMrDE4Z8KgTsS3DHPbrzq/XKnpefxxtMDnVqCudIwJd4t3rws3qKt5yfQBisvUI2XlZFVJC1HKvy0vMsH6JHuku5u0BPadLOFAdom+R1cEmkPqgr4ni1nX9OHpZWlu2o1E3/f7g5HA06XW/dur9zWqeh4M2qkfAWXcc/LllgyKtiyn4NfP3NrV8BhkpTNP7/ho4wvrNCY5v2CumXucCKBdixASWob6PQ4FiAqt39x7o/oDL5vfRxFVQumQggdsdzHHfq1b3xhRzALLC5wlfJcJ0yQ4NEnbqCrHILtD8soPexq0pVFVhaOVHnZO/D3Pb64JPc8UAlLws4oqqDi2x6s+XdGWZZTt17cfeaXV59xzfnmdfzyBvjXILbWRke8dK98i1OQV8KAKeG5sX47ml7m1t4600LVNRI3jv11dzrPmeW5tNqlyr+1evtIuYI59HO/WUKfPF0IRHC9wH3tuSd3zDjpJaBmKw6HxmYeHnfPU4NyurZk5urfP1dp5lBe01jcdASn4fxzJ91LznTiP82oyFDkNTxOGdMKs6iItEYbV2xVVgXJbVZLLtlEWBiXFMGGwrnZ7CluX1hHccn6Xao4FJfepX/KEeaFba4U0caftfg62HVnZIniphjp98Xj7IMRFhTTK1sVJRn45DqlrE0633ye++oNpn2zm6WVp/LIvm1nLd/gU6kBKc326CEjBP1FU7nWMB4CiYoqOo0Woya+glA4tw3Bo+p5Yq/zfFQH0iY92O/7LKrKy6VAeCzYc5qmv07yO/yTw9Dc7jHYFGNUrDkVIHjF9zgPmxW79y6SF22wP8aOWzDmJMS4fqF6nb6H9Irf+CSKbhZbnaYu7wB0vKPcdh+mDuqX51JOPfL75CAs2HOaHnZlYHbJGoXbVqJzVjIM0LgEn+As3HEZKaH31I+4XFBMh7bsjFJUurWtWtZ0YX1ChC6znimk2KYwf1Am1Gq8ehya9hMhq17C52A8kkFlQxkzzv7nb5B4bVSJDmGp9lJ+1foSYFfrERyPcXkvwpP1WvvSo05ekZLLA8gKxVLnuSunbb8EXg5JiMCkCgbvDkgCj3ZWOMWFGNKRrX889vnNfDxhORM5qxkF1v3EJOMH/6OcDempqzb06bmjH3rT7+8sAnCj0r4KNU+U/r6vvxJzXJScwaUgnZo3pg1rNMukqDiZVMH5QJ0xunTVuyHydqYp7SasiGcYU62NsUXrp6a8r/dk1DxVCQ+Fh2x0sd7jnTummHOXflheJRvdjEMJ9JVeEfqTXs12UW7uqwKNX9GTWmD6VglnVP8SsMGtMH87v5v73cM3G63yfE4d0cjvK89zX7zhWgCbdNYOgwa/xCLzjvMqzLu8Em1VHeZ1a+V+zLjkxhvsvPZvfXGL8Qf9yO4NbJg3pRPd2USzZksHGg7nsyyr2eo4AerdvYfwsAAUHL5vnMk79ya1vnozkRusMdtIFRUq2ZxTwx9E/uKZ/PGaT4nXiUFWnz8ZIdYvR3lM5zOLI/+Ha4kcpklXvWREw7YIuzLiyJ6BrSU9Veuipir5W7DhW4HaO37dDVe797u2ijNRdEm9vxtvO62w824nnvl6Cm3OOq4NRMANPwwk4wb/lvM48/tUf2HKqz7P36BU9PW+rkeTEGMandGTBhsOALrTjUzq6fTFdo9Umvr/eSzglsD2jgB3HKgUMO3PM/8todb1bv5OyBTdaH2d35Zm8q9vt0q3H+NuAeNKOFXpNLnZMTLfdx/u8ynC1Ki12N/tePra8zBTrDEoJNZ754S8HOZBdQuuoELKLKgwht9o1bv54IxaTu7LYtkUo6w/ksOdEEXmlVkMDKSqz8e4692y+H/x8kJG927n9fTy98MYN1MubObP45pVag+65jUjACf6kIZ14/Ks/agzOqQ/XDkxwCxe9tpoElcmJMTx9dW83d1cn+uooCVPszFH/xSiXKroAWbRiovUx9svqcwZ8vfVYtZ6EFViYZnuAeeJlhiq7jPYU5U8+ML/CzS51+mwOaeQk9KSw3O7VtmpXJit3ZlaV6FIFi6YNA2DL4Tw3px9fPvvOYz6noINeuttq19h0KJeZo3sH3XMbkYATfKimco6Lql+f1aQuOd/zSq1eQm9S9FU70mTny1b/y1kF7kKfIVsz2foE6bLm3Pm12efKCeFW60N8YnmJZGWv0X6uupP3eJ1ptgfqVafP8/1YHZJ31+7np70nvbQbk4+8fqnpeW6CPm5ggtsKX58EoEGqJ+CMews3HK4xOEdQ8xl+TXh6+FXH0C6xXuf0AzvF8OglHfmpwzucVeCu3h/S4rihYqYh9CEelsJIi0qrcDNto/w79iohjJutj/CHluTWfqG6jbfMb2DCe0WvD6npuVTYNK9JwZcNxXOPn1VUgSIECnXL2BPEPwJO8D/65SC2nCNuba6VcwYl1c/329PiXJMFOjkxxuvIULUWMS39IVqccE9ZuF/GM946k2NUWcorHO7HgMVWB7mlNrKK/Pe0KySCG62PsUtzrxE4Sk3ldfP/olQT9FMXcktsPjWQfVnFTJj7G49/9Qep6XleMf2qIlj750k0KVEUwdRhSXy5JYOJc38LevM1Eo1SQquyzw3A0+ja5jYp5aRGHGejcbKo3CvdlquaX1fDHnhni5k5urehtlZngb7l/C5G7bloinnL8S844l68MjuiG8VXL6L/xgKv/XZcixBOFPpXnDMqRKXIRwRgPlHcaH2czyzP0lWpSpx8tboeK2Yest2BbKK1weaQLNxwmM82HjYyAplUPaZfUBUpiJR88PNBw1EKgsa9xqDWT9WlhNYVQC9gohCil0efbsBjwHlSyt7A/U0w1kahwqb5qJyjC77qkWrKXzzV1O/SjtfqcjppSCfuHN6FAa2srIh5mdgCd6En/hxa3/MD9rDW5PuIG2jXIrRa3wBPfAm9k2yimWx9nHTNPe3XOPUnnjN9TNNmHdRtA7ZKTz67Q6NDyzC3LL1OJx7nKHw5/gSpO/5M50YJLSmlFXCW0HLlduBtKWUegJQyizMUuya9VH1nZl2ljjXnnQztEotJVfRgHVXhij7ta41AS03P49vftvA/xU/Qtmyf+8WOQ2DK16SeFIx/71efYbBbMwporFwgmbRikvUJjmqt3Nonm1Yx0/QpTS38TjQJMeEWtyy9s8b0qcqya1K8HH+C1I/GKqF1NoAQ4hf07cDTUsrvPR90JpTQ0qT0kVlXX/Ftmqx/kQZZdaDevV1UrRbotJ1pfKI8Q5JntFzSBTDxMwiJZP2BfdgbvtX2i/0ZJ7nk51J+HKcRH1W1Htxi+p5yLLxsH0/dvfTrhgKkHSswkmc68+p1bxfFl1v0fH+uqb2C1B9/BN+fElomoBtwIXqlnZ+EEH2klPluN0k5F5gLegaeOo+2EXDYa66cU5+94/oDOdgr96AOTXoln3Dimi3mmuNvEuMh9AUdRrAo4QVabM0lr/QERWXuSTObAumwkf/LIgrXLyZTalzyicLaqeG0jagS/rtNyyiTFt50XFvDk+qOQM9k/GdWkZ7dR4HFqRnYHXr6r+tTOhr+EEu2ZGC1ayxOzeC65ASvCSA1PY8lWzLILqqgdVRIcIKohcYqoZUBrJdS2oCDQog96BPBpkYZZSORmp7nXTknqrVROae+R3n+5H5LTdez1docErMq+GLKS4QWHiIsfw8A+Z1GcsHBKRTvTzfKUzfxAos1+zA5y1/FmrnfaNudrTHy01JWT42klUtOkQfNiynHwvuO0XV+Hdd03qqAdi3DOJ5fhiZxyz7k0MCh6Z+NtdL4t2RLBhd0a2Nk67XaNRZtOMyXWzLc8vpPnPubW9KQxZuPsGjasKDwV0NjldBaClwEIIRoja76H+AM4721+7HleuzvXVb7fgnR9fqiOPek4wd38iop7TzWe3ftfqyVaa2sDskXO0sIu205xHaF3mNZlPgsJXbVUKV8Rfs1FlJqFG7+muPz/uEm9E4OJlzFFMdTFEr3UmJPmBcyRf1vnV6rZ7sonh3Th65tI2kdaUECx/LK/MpGLNGNsT/uznJTMSVQbtOY9c0OQ4vyzHxsc8hgYE8N1LriSyntQghnCS0V+MhZQgvYLKVcVnltlBBiJ+AAHpZSnnHZEw6cLK52fw8wflDD7A5fVqqjSypXI8Crco+T7KIKiGwLt/yX1CzJka0nMKmKEaDiCyk1irf9QESvC1Es/qf4csVeeJKcb+dQnr7N65oa2YrYK+8nrPNA0oCp1kf51PIiEaLq2HCWeT4VWPjccZHX/b7YfaLIOLasD87tkxPXXMLbMgqYMPc3LuzeFlUVbkFSZlUEA3tqoFFKaFXWy3ug8t8ZS16J1Ts4J7Yqj76zqER98JU15mh+maGietKmMr9darbK5I/0L2dNhwr2wiyy/zOHisPbsWYdJHbUXXUeY8nONeT88A6yosTrWniPC2g16m7UsKq/wRZ5NrfZHuJj88tuRTpeNH1AhTSzVDvf6zmeNLbS4vk8m0OycmcmZlUwqlec0UcAn2867LPyT5AA8tVPTc8j0Z4JVQAAIABJREFUp9RWraovqZ9hz4nnPj8m3MK/Vu31+cVXFYGkythnlPXy6Ny1TQQIwY7fN5P5f//PENji3/9D+FmDCDsrxa/6fo6yInJXvEPprnVe10RIBLGj7iK85wg9T4EHv2m9ucP2AHOrqdP3neZ5wHPqkegTwOb0XBQhyCm2ev1NnJV/gugEjMuusxikl/OOc8UX9ffRB9zOnhfcNpS8UqvPHHwASGkkoYwJt2AxKT7teEO6xNK5dQTmNomo4e5lr3K++xeOskIfd7lTdvB3jn803afQhyb2I/6Wt4jodaFPoXeyVuvPdNt92GRVQU9VSN4wv8XFypZq7zsVOJ2YJLqLcLYPoQeoUZ0KQAJG8IvKbGilBd7BOVG6D/wdF3RpsBroGqTjmkbaM/WWMwmlza6RdqyAvh2iiQxVvYR/+fZjrNqViWIJo/VVD+rnXc5nlOSR+/1baNVYADVbBbkr3yPri6dwFHuYW1QzMRffTtvxz2Fq0cav97ZCS+F+jzp9ZuHgHfMczlfqv4dvKP46MTkcwaSdrgSM4P92IAdHaQGm6Dic52SuwTmeGWEagmt55wmDO3mZ542inQIWbTjMpkN5FJU7vFaqgjK7Yf0O6dCD6GE3uF0v/fNXSnas9nr9ihP7OD7vHxSlfuN1zdy2C+1vmkOLQWOM9+4v/9GG8rDtDjQX4Q8Rdt43v8oQsauGO+tG/k8LyP52DiU71+Aoqb81vme7KEyKM/15UNV3JWD2+Fa7hqVNIh3u/BDNVoE9/zjSqufWM9WlxG0tuJ4pK0CXtpHeq5IQSCmpbidQHdHnTqDswGasJ6pcfHNXvEtoxz6YotsiNQeF6xeT/8tC0Dz98wUtho6j5fmTEWrd4+2dfKkNJ8Ru40Xzh0ZbmLDyoeV/mGKdwRZ5dr2f7aRk1xrseccp+WMlAHGTXiK0Y586PyehVTh/ZhUjkTgk7DlRFDTuVRIwK77N5UhIMYdgaZNESIceALSMqL8guJKansesb3YYjiQa+Myv55kQ01+EaqL16AcRpqq4e2ktJfs/r2HNPUrmgkfJ/+lTL6FXo+OIm/wSMSOmNkjonSxyXMIzthvd2iJFOfMsL9NbHGzQs+0Fmdjzjlc1KCYscV3r9awVOzONo0CHJnmqDmW6oXlX8wkYwS+uwf21c6x/6bRrwhmauz2joMZ+qsAji27dMMd2pOWFN7u1VRxJ4/iH91BxbLdX/8h+o4i/+U1CE3rX+zWduCpGHzuuYLZtgtt1Z52+s8UR6kvZoa1uv4d06FFvnwVPHJrkyy3e9RR8CXhzr+YTEIK/cMNhThRVH7t+dlz9z++dOI/lalvLnWGoDSFq4FWEJp3j8WD3rDlKeDRtrn2S2CvuQwnxP2twTXgqKu84ruFfPuv0PU9C9m8Ubf0em0dcRG2Uewh+aNKAeo21OtbsyTISgED1At7cq/kExB7/u7TjNV7vHR/t97OqK8vseo6vKoIBHVv6DKdtDIcWIRTCe15QKSTeTww7axCxV9yHGtH0+9nX7eMIwcqdpuVGWxtRyOex73HBkhwO/ldiatme0KT+hCYOIDSxH2pYC5/PklLz8igMS2xcwT+aX87CDYdZvPkIT1/Th+/SjhtOPhU23evS9VSmuSb3DAjBv6JPe37a610vzklNBTJd8cy04+oC6plsE+CG936rtYZ8XdEqSsld9T4lf6zweT2080DajJtZ47l84yJ4yT6RUKxMNf1gtHYIs/LjTREM/7iEI/nHKd56nOKt3wMCS7uz9EkgaQChCb0Mm4Ut6xCai2+CCInA0r5bk4za5pDM/DrNLcmHRI8O7BMf7ZYivDkm9wwIwZ80pBNPfPVHtautv7O5L/XPM0W0s9+2I/mNLvTlGTvIWf4a9gLfaa8Byg/+TkXGjnpZweuPXqcvBBsTTVXHi0ktFVZNCWfEvFKOF1eJl/XEPqwn9lG4YTHCZCGkQy9CkwbgKM51e2pop74IRaUpkOhJWfTRV+lNdrtmTAjN2b8/IAR/4YbD1Qp9q3Cz3x+smzrvUvvdswyU1a41akStdNjI/3kBheuX4KXaCwVhDkVanWWqJNnLXyP+lrcabW/v1xhReMJ+K6HCylj1F6O9W6zKykrhzy71/hSk3Up5+lbK07d6XbO07dKkY3aiKqAoCg6HBkIYE0Jz9u9v9sa91PQ8Zn6dVu31/DKb3xZbozx2pVOO0+124YbDvL16n5EsQpP+e5TVhvVkOsc/eZDC9YvxFHpTqwTaTXmNNtc+4dbuKMwid+XcxhlAHdBQeMh2J/9xDHZr79VGZdUd7WnbMbFOzyv4ZSFH372VnO/fpGTXTzhKaz4xqS8X94hj0e365ypc/saqj/z/zYVmv+I7s+NUh5R1C85JTowxnqlJ3SD01FJ9G6EqotFi6KXUKNq8jLy188HhfRQZlXw1LUdMRTHrUX5RKWMo2vy1cb0kbSXhXQcT3v3cxhmQn+h1+qZjYY5bnb5+LYr579QuTMp/kpOH/qTs0FbKD23FUXSyxufZCzIp3vZfirfpeQAscWcRmtif0KRzCEnoZbz/hrBqdyato0IQVJ1cCPSip81xtYcAEPzaZuz6WGyHdonFpAgjsYZzddcaaZnXY+Zfpzx9u9c1NTKW2Kv+SZjHMVfMiJsoP/g7tpzDRlvOf98ipENP1Ej/v7yicsPbkHdiq6ZO3wDlAPNbvsdNvWZQ2muEHjSVd4zy9G0Ubv7aK4DKF9bM/Vgz91O48UtQzYQm9DQMhZa4s+plE3Bo+nbQogpMqq7yq6quDNc7B+MZTrNW9Z1Hb+YaHGaerqzwWheSE2O4PqVj7R3riJSSkp1rOPbRdJ9CH95zOO1vfdtL6AGEyULrqx8EpWou18oKyfnuX/+/vXMPj6o8F/3vnclMIhDImAsJhiSCgCAUJcjNC17ACiqIuKuA9VB2q7CtPVbbLcVujgcv29O9u2trPSja7l37bFBbFalHuxW3Qr2EAsolgXILBMNNAjEgECaZ+c4fa9ZkrZk1l2Qmk8ywfs/D8zCTNTPfN7Pe73u/94qKUw0RYGBBTyYP60tx78R2Ur1P3zr/xabnLw/06cvGi4jgOv8Cci+biiPbHESVld8fnDH2JV8LzXVb+GrtSxx+6UEOLP0OKixUOX68PsXEwYWmo9ysZZ/yiMHvnylkrOAbAzOiBczE68oL5bZRpUmN8fedOUnDqp/R8Kd/DSuU4cjuScEtP6Zw2j/izOkV8T3cfQeSd9Uc03Nnajfw9aZ34hqDAnYfPcW7245wypt4G61mspnn/TGf+c0htxOc23je9QvcaEcYX/PXpvwDgOI7n6T//3yZom8tofeY23D3HRjz81z5/RP2AqzZeTRo8ferttp/mRa9l7GCbypwEYFEMrYqyz0smT48KcKv5czfx+m//SXsbznlIymZ92t6DpsY13v1HnMb2RcY+p04XabiovEiIX6Jjk7zFOcx1/uwZZ++Z1zPkEUrZ/dvMY3RVVCOs5cHhyuH8y4chefaeZTM/SWl9/8nBdMeptfIb+LsE948NBlRfq0+PzUHmshytH0Dxgo+Oukex5+UFloiMhf4F0A/pP1aKfViEscZFatoOt311twS+aa//uKiDp/fNtY10njay5Lpw6k52MSKdfvb3W3O39LMV2t+Z5k+i9OF55q55Fbe0q70WXE4yb/5QQ79+/1k5RVTcPNDuAsrLK91SHgYrs6J5lYcAkW52dx66QVMvqSYqtpjnDzTwrK/1MZVLDP4XoE+fS+7H+diR1sc/zedG/g3lvLtfeZdOpIAO3v0oefQq+g59CrNPvDVYc0VuG8TzXWbyUkwyk837m2pb8KV5WDSsCLW7DyKz2eO3osWyJUuxBR8QwutyWhltNeLyCql1LaQS19RSn2/E8YYlUg/gu56+/vfreer09YJOvdOtFYf9YXE08NtGbll1SvP4ZR2GffOHtpFw1s/p/V4eNKIu+9A8m9+CHdBx4p/uvKK6XvnE7gLL0SyrLPxnA7hlm+U8ObmgxE9EX4F1w/ty8KpQ1m+bj/v1hxm64Gmdgm9zlfkcpd3Ea+4lzDQ0RZCPc35Kc8MymLu520GxXh2bhHB5SnB5Skh99Ip2tk+gWhFl1MYVtKbLfVNmsHW5+fS/nnMnzgwbFOJFciVDsSz4wdbaAGIiN5CK1Twu4RoP0JluYeK83uwycL/e+ul/Sx/LKNQ+5W2C2S7NOGuOajdFLuOnAxqEi2tfn77Ua2pwms0lN9HU9UfaPp4RXjOvDjoPe528q6YlXD6bHZJ5Lz4i4p6Me+KC1n8ZnVM92P1gSYeePlzVm4KbaXQfhrow2zvI7zqXkK5o63L2t1DWzl1Uw7/8P+aweHsUNRhomf7/z1tOEOKc4O9D6CtnVfofZIJcfzJaqEFMFNErgZ2Aj9USoXlZnZGC61YP8LhE81hrxlanMvTd14W9jyE2wb0Jg7/9Ga1ZQiuAvYcDa9aa4Xv60a+fONxvAd3hP0tK6+Y/JseJKd0mMUrk8uki4toPO21LNvlcmpFQvSK4Fvqm9gcI9W4PRzhfOa0PMIr7iVcIG1n5gWj3TS3Kn5SMzDY4CSV1BxsYkhxLr6AK9On4NFV1Qwpzg0T/NC8jHTb7SE+4148LbT+BFQopb4BrAZ+Z/VGSqllSqnRSqnRhYXx1XqLRWiRy9AfwWrw0cpojxuQH1YjT0Qixt372+HzduT0QrWEpwf3+sYNlMz9VUqEHuDFj/ZaFvkU4O9G9+eVeydw1aACHNI57TLrVSFzvIs4ovJMz/9wXDZPXJ9Dqpp0GlFofRGMv7PelMMKY33FdCQewY/ZQkspdUwppd/RLwCVyRlefET7Ec5aNLN4c9NBlq/bH/a8jrFCjgDXXVwUthh0BMlyUXDLj4L+aS1n/p+SmjMfD36lgtln3yjtg8spOEU70ug95x6YNBh3liN4gwjaEdpBcjp77VMlzPEu4pgyL8LOwthuu/ZSfn6PqF4Jp8Dwfn3ClhtHgpWXuzPxqPrBFlpoVvs7gdnGC0SkRCmlW2ymAcmrvJgAG+saaWoO90crYPGb1mpcVe0xk/HK6RDmTxzItUOK+OnKrR0ybBlxF1bgufpumr+oIf/G+3H2zIv9oiTiCLSb9vRws+StGrytWoPKO8eUcZuh0aRRnTUaOQGeXr2Tj3c3JPxd7Fal3OVdxAr34+TJKX7Wcgf/1xnagT1x6o5rCUwCjK7wsLm+CZ/Pj8Mh+AP++iVv1bD45ktwO4UWn8Ih8NitI9J2R49Fslpo/UBEpgGtwHFgbieOOW6qao9FVtEDXW2tDDfZLgfeFu3GWDJ9ePCaLId2UwiEue4cAiV9cmj42mupZRjJvfxWci+fkcKceQ0B7hxTxvB+fXin+lBbQpFf0S/vPMuzrNWN/8Ckwazbe9zUGszp0PIe2rsYbFfl3O1dyEjHHn7vu6Ej04obBQzqm8vCKUODnY5e/uv+YBGO6oNNrLhnfFqf3eMlWS20fgL8JLlDS5xxA/IjdprJimCNtSqo8ewHuzn41ZlgK2wRyAqc+/X39iutuks8tLesdSycDijunRPX5zecPMviN6tNueh65x+9L32sG76y3MPtlaWsCKQ7C5qqvPVAx4yAW9RAtviSr+JbIbTVTXj9s/pgLIMCXt3wBTNHlXLftR0r7plOZHSSTmW5h4uLc02tmINE8WPpO93ydftZHLDmu0ISOK4ZXIhCC/G0aoqZSvx+GNavj0nwx1R42FTfZBqbwyGs3nbEpK2IwNzxFUG1P56AlI11jdqC4dQWP1eWg/ED8mNa/11OIb+nm8MnItc/7EwcopVZM/6uRqWr1aeCpbfAHBgGZJQmkNGCD9D//B6Wgt/qUzy9eicPTBoc0Z9v3Bm9PsXkYUUI8P72I7y37QjZLgeP3nIJv3p/Z5fdzKDtVk2nvUHtRoCJQ4p4OKDSnjzTwosf7Y2Ynlxz6ETcASnGOIcsp4M7xvRn5qhSqmqPxezj1+pT9O2dY/lduZ1i6m+vE09vwHhf41eweFU1yq+CGZWh67++Dpjm6RCtQIcvfSP1QsnYWP1Y+IGPdzdETL7QjHzmu+K//3aE97YfCfp6vYFzYdn5qbPIR+L4qbaecQqtZZju7cg9z2Xpsxc0Q9+U4SXBdl+xAlKMcQ4+n58LArYB3TbiFMhxOZh/9QBGlvbBmBjpdAp3XF6G2yJb0krooWPBeLEWH5PxVjRNRP8ubhtVGjbPFp+iJcMq7mb8jl+Qa04vHVqcS0FudtAqHWmH0wODjG2uQzvfKOAPG75IuFx2LByi7Uwq5DnjDTygsBd7jp4KXvPCR3spy+9J42lv0Gevlwy7vbI0WFBSV12HFOfGpcpGCpiKFNSyfN3+YPCTQ4QhxbmsuGc833tpPcdPRe51oJPksoVkOQUHWvadQzTjrdXcQ6smIxIWs5/OZLzgD++n7To+pamTj88YAcD6fcejhlzqN/Lrn9Xzhw1fWO5IisRr5Meid04W//6dMcFxaOdSMe3gNwzry70TB/L+374M6xyjAkUjY1WMjWTBt7ouUtSa1Xs0nvYG6wG0tGrlq5+cMYKLCnvx1wT64sWLHp6sh1vPNOzoxvHHis6zek06k9GCrxtx/EpzxT06rc01F0/IZWW5hx2HT1Jd0kTDKS8HGs8kbWzuLAfzJlTw3NraqNddF8ggrCz3cFvgLK374PWF696JA6ks93DdxUW8t62tAq++CLS0+mk87U2atTreRQLCqxXp5as3ffFVuz9XAE9PV0RNYUyFh7OtfpORccyF5zN7bHh4eLyLnFUV5UwgYwU/1DinR6vpxHPzPvX2dpNgdsTQFBGlOHk2PLjI7XTgNZwp3q4+zLcD5Z+MY9bVU08Pd/DMOX/iQNbs+JIWnzJVju1K9VSvVqRXOvb5/LxTfcjkTpw0rC9FudkcPXmWDwPjD/2eBZg1VotBWPSGdVvus63+MO/C8HY0SzmXyFjBDw3ecUj7im5srGtk2V/Mu3EylXqfX/GlRVsvb4ghIZINQn8cmh78d6P7R1Vpu4LbRpXy2mf1QQ1lyvAS01FrfkBj2VjXSKFhAWj1Kfy0RRvqc3IZouuMto+tB5qoOXgiuEA7pOMVljKdjBH80GIcupX5bCB9dlRZeGhspHZYEB66G4tY2oAAk4f1Dd7QAB/s+DLKKzQcUaoEGS3PZ1v8/HTlVpQyx9x3B/XUyi4QalALrXHw6LThQcOk0Tbx7Ae72xZ0pX0/vqBWpz2pVTtWMTWdaL9/ppMRgh+pGMfimy8JWpT/uk/rW7/invHBG23WC1XBXWfF98y+2XED8qNWqDH+ze0Upo4oiZqz7nBAYW428664sM2nHsMw6HQIj03XctOf/WB3mBDolmdvix8/bT5pb0v3Kw5hdV42Pg6tqxDJJmG0tocaOXWXnJUhM1TIM6GKTiJkhOBHKsZRc7DJMs2ystzD64HmF6Dl2z+/Zg8j++cFb4zKcg/3XDXA0vimC6TRUhzLt6uXcHY6wlN8HYTH/jsFHgu4mqwKg+g36n9+dxxPr97JR7saghpHNC0hlO6y68Vb3CI0eUg3cupuSl3T0eelEyrkmVBFJxEyQvCtbprl6/bz8l/NqbcuZ5tAhO61720/wurtR0yr/8KpQwF4fm2t6Xrl1wyFTwRcgwA7Dp+0FGojCiyj5yYN68vaXUfD6gM2nvZaFwZp8ZuiDh+YNJj1+7SkGd03Hc9N3J12vfYUt7Aycobu7sZ5zRxVGibk4wbkk+V0aM+jBTydS2SE4Fv5XBe/WW3SpEeW9mGxoYb+zFGl/DEQfKNoMxJ5Q1b/hVOHMvmSYp5fs4f3tx9BKXC7zDvSxrpGlrxVg9+vnS+Le2fHlTCj7973ThzIvRMHmnz1xl3PneUw7fh61OH6fcdNO397d+6u2vUiaRkdsUlYvSZ0XgostQmfX1tolSKo2emLfaaTEYIP5htg0RtbTTtrlkNMQq9fv+Ke8Ty9eqephbaV9b+y3MO9EwdSmJsdVO2N7/XaZ/XBCD+lFIeb4hP6KwcVmHIFjL56o1AYVdt3qg9ZRh12RGi6onZcKrSM0HnNHFUaPI6ZjYTm1z23tpYV6/dT2CubeVcOsPT/ZwoZI/g6G+sa+ePGtsq1TiGi6mtSkw3594ApRdVKddQ/67XALq0vM8b+a9FwOsUyQchKgPXHVbXHwlxhiQhrRzWFREiFlhFpXqHG2yynhBVJbTrTStOZVha9sZXn1+zh3+64NCPP/hkn+FW1x2gNLOV64YloK7du/X+n+hBThpdolVZDrP2hbrPXPtMWljkvVpli+QWtHPXaXUeDMd5WYb0CXDekKGh8Mt5YVmqwVTnvaOG37SHVLr9UaRmx5lVZ7uGVe8Yza9mnEROE6o6fZubST7j10n4Ri7OmK2kt+NEaaeg31m2jStlY18jrn9Vbqun6+fxsi59P9xxjVFmeydr/+mf1wXZZxrBT/e9GoTee1/VxvVdzmJWbDuDp4WZPw6lgPv+HO74MMybqAn62RVs0lkwfzuyxZXG7uqJ9L92FrtAyoo1l3hUXxgybXrnpICs3HcyoBSBtBT9WIw2joW/WC1VBYf7jhi+CvnzQNAR91271K9aHpOgqtBvkmiFFvBuIg/f5/Ow+ctKk0k8OJMoY1crl6/YHb6rDJ84y/+oB5J7nCpZ8ClV3Q8ei1wVszy7ZnSz1kegugUWgGfOqao+xKY4S4is3HeTdmsP8vht+p+0lrnx8EblRRHaIyG4RWRjluttFRInI6OQN0Rqrs6JOZbknuCMu+VONqQpNaMnksHLaSrMLCFpgzsyAxvChMcpOhPX72hYIAUb2zwu7Gd6pPmR6XHPoBPddexEzR5Va5r+HjsWv2uIOopUQj/d7sbFm5fev5NZL+8W8rnHtSxzbsZ6ZSz+JWqU5HUhaCy0RyQV+AKzrjIGGEmsX3FjXaHl+M/ryoa35pZbFZ53C+uwHu01eAr/fnERiVYZ5Y10j57nM3V2mDC8JfmYk41PoWIz57vHsMpnQ5aUrePrOy/j2+AqeX7MnqNkZObP3c058+iongJ7Dr+fhM9/l8bdq0nb3T2YLrceAnwE/SuoIIxDrrFhVeyzMqBbqy9eZPbYsaiEKTw+3Sa0PNQU5Qoq2G9VtvSfbHZebjYyRBDnWWGLRnc7Q6UZluYdld49m+br9pgxAv/cMx/78TPDxqer3ObN3I+dPXsDMpX4uLe3Dyu9f2RVD7jBJaaElIpcB/ZVSb4lIRMFPdgutaLvguAH5uAx13NxZDkuhj+e9qg+2nf8EghF6+gIQWqrbqG6LX3HDJcXt8gknegbuTmfodET/rR55YysKaK7bjO9kg+ka/6mvaFj5z5wePIGNkxdQsbCJJ2eMSBvff8IttESrFf0L4KFYb9QZLbQioQfozB5bxpyxZWFJOPESGhfgcmrW9lljyyLWqdPVba1ik+Dp4U7KnGxSx+yxZex96iauHlRAj0HjKL7rX3AVlIddd3rnJxz8zQJObvovHvj5fzBs8Z8tazh2N0TFaJcqIuOBR5VS3ww8/gmAUuqfA4/7AHuArwMvKUZrqjFNKbUh0vuOHj1abdgQ8c/dhmc/2M3P390RDJedPbYsGKMfzW1mLOFsTKqxST901V/5Wmj69FWaPv0D+MOLqACcd9E48qfcz/+Zc2W32P1FZKNSKszYnnALLaVUE1Bg+KAPgR9FE/p0wiouQCeaSq13o1Wcm9lfmcTssZp9Zvm6/Sxyuugx5AqOvfNLvId2hV17ZncVB57bxIONj/HoqmHcPrp/WOxIdyCmqq+UagX0FlrbgVf1FlqBtlkZTXtcaUb0BSOektU26cHssWXse+om3IUVFN/1r+RdMw+crrDrVEszjf/9IqeO1rN83X6+9fwn3U79j6nqdxbpouonQneOoLNJjMGPvI3XpzixcRWNq5dZX+R0kXflHHqPmUGvHDc1S25M7SCJrOrbgm9j00Geens7S9fs4fTf/sKxPz+D8lpXYXYXX0T+lB/gLhqQcsu/Lfg2Np3E5Y+/x+GjDTS+/wK+r4/TXLc5/CKHkz5jb6fPhDvpX5DLRwuvT8nYbMG3selEgpZ/pTi17UMaVy/D3xzes9GVX0bB9H/EXVjBawsmdPoRMJLgn7O982xskolu+Btc1Itel1xLv+8upcfFV4Vd13qyAUdOLwBmLv2Eu3+Tkgj3MGzBt7FJIu89dA2vLZhAVs88Cqc/TOGMR3D2bNvVPdfOIys36P1m7a4GJv/8w5SP0xZ8G5skU1nuYe9TNzGosCc9Bo+n5LtL6TliEjnlI+k18pth1+86eoohj7ydUpeffca3selEnnp7O7/5eK9W1LW1BckK9/sbSXaxD/uMb2PTBSycOpRdT0xl/tUDYgo9aMU+nnp7e6ePyxZ8G5sUsHDqUC3qz2mV82bm+bW1TP/1R51a7MMWfBubFLLziak8OWME2VmRRU8Bm+ubWPTGVqY8vbZTzv624NvYpJjZY8vY8fgU5l89IOa12w+f5FvPfcIjb2xN6gJgC76NTRexcOpQXlswgdljyyyLXuj4lBYgNOfFqqQJvy34NjZdSGW5hydnjOCPCyYwpiJyFJ+xZ2IyhN8WfBubbkBluYdX50/gyRkjuKiwJ0W5bi7Iy+HyCo9WzYm2nonJ2Pm7tR/fTmu1sQk0fflTDZsDtf8dAldcVGDZgi2URCrwdAnp0BjCxiZVbD90Ivh/v4KPdpm7JbeXbqvq240hbGw0qmqPmfo6AKaSbh2h2wq+XbrKxkbDKAvuLAdupyQsF3Gd8UXkRuCXgBN4USn1VMjf5wP3AT60arv3hHbaCcU+49vYxI9RFoC45aLDhTgCLbR2YmihBcwyCraI9FZKnQj8fxrwD0qpqAXG7CQdG5vOJ5EknWALLaWUF9BbaAXRhT5AT8K7TNnY2HQjktJCC0BE7gMeBNzAdVbo4b/AAAADVElEQVRvlOwWWjY2Nh0j4RZawSeUelYpNRB4GPip1RulsoWWjY1NZOIR/Hqgv+FxKXAwyvUvA7cmMigbG5vOJR7BD7bQEhE3WgutVcYLRGSQ4eFNQHhvIRsbm25DzDO+UqpVRPQWWk7gt3oLLWCDUmoV8H0RmQS0AI3A/+jMQdvY2CRGl8Xqi8hRoC7FH1sANMS8Kv2w55VepHJe5UqpMINalwl+VyAiG6x8mumOPa/0ojvMq9uG7NrY2HQetuDb2JyDnGuCH6Gfcdpjzyu96PJ5nVNnfBsbG41zbce3sbHBFnwbm3OSjBR8EblRRHaIyG4RWRjluttFRIlIWriM4pmXiHxLRLaJSI2ILE/1GDtCrHmJSJmIfCAin4vIFhGZ2hXjbC8i8lsR+VJEqiP8XUTkV4F5bxGRUSkbnFIqo/6hRRfuAQagZQpuBoZZXJcLrAWqgNFdPe5kzAsYBHwOeAKPi7p63Ema1zJgQeD/w4B9XT3uOOd2NTAKqI7w96nAO2iJcOOAdakaWybu+DHrBwR4DPgZ0JzKwSVAPPP6HvCsUqoRQCn1ZYrH2BHimZcCegf+34foSWLdBqXUWuB4lEumAy8pjSogT0RKUjG2TBR8q/oBFxgvEJHLgP5KqbdSObAEiTkvYDAwWEQ+FpGqQMm07k4883oUuEtE6oG3gftTM7ROJ565dwqZKPhR6weIiAP4BfBQykaUHOKpi5CFpu5fA8wCXhSRvE4eV6LEM69ZwH8opUrR1OPfB37HdCeuWhedQSZ8eaHEqh+QCwwHPhSRfWhnq1VpYOCLpy5CPfCmUqpFKbUX2IG2EHRn4pnX3wOvAiilPgVy0BJd0p321rpIGpko+FHrByilmpRSBUqpCqVUBZpxb5pSqrtX/oxZFwFYCVwLICIFaKp/bUpH2X7imdd+4HoAERmKJvhHUzrKzmEVcHfAuj8OaFJKHUrFB3fbTjodRcVXPyDtiHNe/wXcICLb0Eqd/1gp1a07kcQ5r4eAF0Tkh2iq8FwVMIt3Z0RkBdqxqyBgn/hfgAtAKfUcmr1iKrAbOA18J2VjS4Pvz8bGJslkoqpvY2MTA1vwbWzOQWzBt7E5B7EF38bmHMQWfBubcxBb8G1szkFswbexOQf5/96gBJ0Dc7liAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "(handle0,) = pl.plot(x[:, 0], x[:, 1], \".\", zorder=-1, label=\"empirical\")\n",
    "\n",
    "for c in range(num_clusters):\n",
    "    (handle1,) = plot_cov(\n",
    "        gnd_mu_[c],\n",
    "        gnd_cov_L_[c] @ gnd_cov_L_[c].T,\n",
    "        color=\"k\",\n",
    "        lw=4,\n",
    "        label=\"gnd\",\n",
    "    )\n",
    "\n",
    "(handle2,) = plot_cov(g_mean, g_cov, color=\"C1\", lw=4, label=\"Gaussian\")\n",
    "\n",
    "pl.legend(\n",
    "    [handle0, handle1, handle2],\n",
    "    [\n",
    "        \"empirical\",\n",
    "        \"gnd loglik={:.3f}\".format(\n",
    "            mixture_full_log_pdf(x, gnd_mu_, np.linalg.inv(gnd_cov_L_)).mean()\n",
    "        ),\n",
    "        \"Gaussian loglik={:.3f}\".format(\n",
    "            Gaussian_log_pdf(\n",
    "                x, g_mean, np.linalg.inv(np.linalg.cholesky(g_cov))\n",
    "            ).mean()\n",
    "        ),\n",
    "    ],\n",
    "    loc=\"upper right\",\n",
    ")\n",
    "pl.axis(\"square\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# full rank"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "1 loglik=0.860 elapsed=0.0s\n",
      "2 loglik=0.860 elapsed=0.0s\n",
      "4 loglik=0.861 elapsed=0.1s\n",
      "8 loglik=0.923 elapsed=0.2s\n",
      "16 loglik=1.662 elapsed=0.5s\n",
      "32 loglik=2.163 elapsed=0.9s\n",
      "64 loglik=2.392 elapsed=1.7s\n",
      "100 loglik=2.392 elapsed=2.6s\n"
     ]
    }
   ],
   "source": [
    "model = GMMModel(x, num_clusters=num_clusters)\n",
    "trainer = GMMTrainer(model)\n",
    "\n",
    "for t, epoch in elapsed(range(100)):\n",
    "    trainer(x)\n",
    "    if (\n",
    "        np.allclose(np.log2(epoch + 1), np.round(np.log2(epoch + 1)))\n",
    "        or epoch + 1 == 100\n",
    "    ):\n",
    "        loglik = model(mx.nd.array(x))[0].mean().asscalar()\n",
    "        print(f\"{epoch+1} loglik={loglik:.3f} elapsed={t:.1f}s\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "mu_ = model.mu_.data().asnumpy()\n",
    "kR_ = model.kR_.data().asnumpy()\n",
    "cov_ = np.linalg.inv(kR_.swapaxes(1, 2) @ kR_)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(0.2941210180312009,\n",
       " 1.0967109330507165,\n",
       " 0.28322551523446604,\n",
       " 1.0858154302539818)"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAP4AAAD4CAYAAADMz1tMAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nOydd3hU1da4332mpCeEAIEQIASQkkAQQhNpKthQLqDSroiI6FWsVxFR+RSvBX+o3KteFRuogN4LgsgHfoqCCEoLUkKvgVASUkmddvbvj5M5mZlMkkkIRTLv8/CQnDnnzJ7MWXuvtfYqQkqJHz9+6hfKpR6AHz9+Lj5+wffjpx7iF3w/fuohfsH346ce4hd8P37qIcZL9caNGjWScXFxl+rt/fipF6SkpGRJKRt7Hr9kgh8XF8fWrVsv1dv78VMvEEKkeTvuV/X9+KmH+AXfj596iF/w/fiph1wyG7++Y7PZSE9Pp7S09FIPxc8VQGBgILGxsZhMJp/O9wv+JSI9PZ2wsDDi4uIQQlzq4fj5EyOlJDs7m/T0dFq3bu3TNX5V/xJRWlpKVFSUX+j9nDdCCKKiomqkPfoFH0hJy+W9NYdIScu9qO/rF3o/dUVNn6V6r+qnpOUy7uONWO0qZqPCgkm96d4q8lIPy4+fC0q9X/E3HsnGaldRJdjsKhuPZF/qIf2pueaaa2r1WlW8+OKLzJ49u7ZD8uOFei/4veOjMBsVDAJMRoXe8VGXekh/an777bcKxxwOR6Wv+bk01HvB794qkgWTejO6Z0u6NI9g5ne7Wbjp+KUellfq2hfx5Zdf0rNnT7p27coDDzyAw+EgNDSUZ555hu7du3PDDTewefNmBg4cSHx8PMuXLwdg3rx5DBs2jJtuuon27dvz0ksv6fcMDQ0FYO3atQwaNIixY8fSuXNnt9cA3njjDTp37kxSUhLTpk0D4KOPPqJHjx4kJSUxcuRIiouL6+Rz+qlIvbfxnfxn6wlsDq0M2Y70XQCM7dXyUg7Jjbr2Rezdu5evv/6aDRs2YDKZeOihh1iwYAFFRUUMHDiQWbNmMXz4cJ5//nl+/PFH9uzZwz333MPtt98OwObNm0lNTSU4OJgePXpw6623kpyc7PYeznM8t5hWrVrFsmXL2LRpE8HBweTk5AAwYsQI7r//fgCef/55PvnkEx555JFaf0Y/leMXfDQ73yn0Tlalnr6sBN+bL+J8BP+nn34iJSWFHj16AFBSUkKTJk0wm83cdNNNAHTu3JmAgABMJhOdO3fm2LFj+vWDBw8mKkozi0aMGMH69esrCH7Pnj297iuvXr2ae++9l+DgYAAaNmwIQGpqKs8//zx5eXkUFhZy44031vrz+akav+Cj2fkmg3AT/oRm4YC20m48kk3v+KhL6u13+iJsdrVOfBFSSu655x5ee+01t+OzZ8/Wt4YURSEgIED/2W636+d5bh95204KCQmp9L29nT9hwgSWLVtGUlIS8+bNY+3atTX6TH58p97b+KDZ+V9N7kPPuEiEAAHM+/0YCzcdZ9zHG3nzh/2M+3jjRd/n9xzjgkm9eXJI+zrZcrz++utZvHgxmZmZAOTk5JCW5jWD0ys//vgjOTk5lJSUsGzZMvr27evztUOGDOHTTz/VbXinql9QUECzZs2w2WwsWLCgBp/GT03xr/hldG8VyYD2TdialouKpk5/veU4FpuKxH2r71JpAN1bRdbZe3bq1Il//OMfDBkyBFVVMZlMvPfeez5ff+2113L33Xdz6NAhxo4dW0HNr4qbbrqJ7du3k5ycjNls5pZbbuHVV1/l5ZdfplevXrRq1YrOnTtTUFBQm4/mxwfEpaqrn5ycLC+nQhwpabks2ZbO4pR0HA4Vg0FBVVXsqva6ySB46fZEXlyeis0hMRkEiyb3qbUg7t27l44dO9bomiKLnSKLnZAAIyEBl27OnjdvHlu3buXdd9+9ZGPwUxFvz5QQIkVKWWFW9q/4uHvMjYpgdM+WZBZYWL0nQz9HVSVr9mdiLfMDWB2SJdvSL9qqX2SxczSrSLePWzcKqSD8l8vE4Ofyp1obXwjxqRAiUwiRWsnrQgjxLyHEISHETiFEt7of5oXF1WPuUCWZBRZ+3peJqy4kJRzNKnK7LqvActHGWGSxI6VEojnHiiz2Cq8fzSoi41wpR7OKKrzuPCfzXKnXa70dr4wJEyb4V/s/Ob449+YBN1Xx+s1Au7J/k4H3z39YF5fIYDOKECiAQRH8vC8Th1ou9gIwmxRaN3L3UjcKC7hoYwwJMCKEQKB50D1X9NxiK6qXicEp1NmFFq8Tgy8Thp8rj2r1QSnlOiFEXBWnDAM+l5qzYKMQooEQopmU8nQdjfGCkpKWy8wVu3GoEiGga4sGbDlW7r1XBIzu2ZKR3WIB+GV/pm7jO4/V5j0dpTZdLXelKnW9QbAJATQINuuvFVns5BVbyS22lZ9YNjG4mgcgkE4dpmxicJ7j1CRcj/u5sqmLb7g5cMLl9/SyY38Kwd94JFv33EsJm4+5b9kJYGS3WN2WXzS5z3l59Z3+hHdvjuZoVpGbrV5ksXPExY6PL3vN074PNBkostgxKILT+aWoHg7aIJOi388p1AKJpi9IN43BqUlQNjkYFH+qcH2gLgTf25PidatACDEZzRygZcvLIyouMtjsfbBlOCTcP38LyXENeWBAm/PeUnP6E1xVcqcQ5hVbce6ySCnJK7Z6XZVP5ZXi/BN7G3ux1cGRrCJiIgLBKdRCEBMRSKnN4XZNSICRZhGBnMorRSI5nV9KoMngX/WvcOoigCcdaOHyeyxwytuJUsq5UspkKWVy48YVavxfEnKLrdWek1Ns44c9Gdz5wW/nncDjjMDzZqt7CrHzd20VFvoMK8uUds/zXWdgKSVZhRU/W26xjZwiK0fOFpFdqDknNX9G+YTjq50/cOBAr70RKjvuC64puDNmzGD16tWA1ochKyurxvd7+umn6dChA126dGH48OHk5eV5PS8uLo7OnTvTtWtXt5iEF154gS5dutC1a1eGDBnCqVPao71v3z769OlDQEDAnzJluC4Efzkwvsy73xvI/zPY985Mt001yL9XJbzwbep5RfA5I/DCg4wVtuSCTAZdeJ0hrafzS/TVGCpRpcrwfM1id7hpEJkFFt0skGiag1PjqMpxeKmYOXMmN9xww3ndY/DgwaSmprJz506uuuqqCiHKrqxZs4bt27e7TVpPP/00O3fuZPv27QwdOpSZM2cCWn7Bv/71L5566qnzGt+lotpvWAixCBgINBJCpAP/A5gApJQfACuBW4BDQDFw74UabF3hum+vViVJXnCokg9/Oczc8b5HqnnSvVUke4vPlNvXF5AdJ8onKZtDdXvtwzlv8H/LFxPXqiUNIqNI6JLE008/za033kCvXr1Ys2YNeXl5fPLJJ/Tr14+SkhLuvfde9uzZQ8eOHSkpKan2/RctWsSrr76KlJJbb72VWbNmAfDJJ58wa9YsYmJiaNeuHQEBARW2CCdMmMDQoUO544479GMlJSUMHz6ckSNH6pl8VTFkyBD95969e7N48eJqr3ElPDxc/7moqEj/vpo0aUKTJk343//93xrd73LBF6/+mGpel8DDdTaii4Drvn1t+GlvBilpuX/aEl0CSN3xB6tXLWfNhk0YheSaXj24uls3fbW32+1s3ryZlStX8tJLL7F69Wref/99goOD2blzJzt37qRbt6pDNk6dOsUzzzxDSkoKkZGRDBkyhGXLltGzZ09efvlltm3bRlhYGNdddx1JSUnVjruwsJDRo0czfvx4xo8fD0C/fv28hvbOnj27grbw6aefMmrUKO9/EyEYMmQIQggeeOABJk+erL/23HPP8fnnnxMREcGaNWuqHeefgctDp7vIOO3sUpta/clecEjOOy32UiGEoFGomT+2bGTgkFvItWr+g2uvu5GcIqtu348YMQKA7t276+m469at49FHHwWgS5cudOnSpcr32rJlCwMHDsTpzxk3bhzr1q0DYMCAAXo67p133smBAweqHfuwYcOYOnUq48aN04/9+uuvPn3uV155BaPR6HatKxs2bCAmJobMzEwGDx5Mhw4d6N+/v37tK6+8wmuvvca7777rVnjkz0q9yM7zVrmmS/OIWt9PCOqkRNelCJaRUqKqsmz7TvvdW0SgMx3XYDBUmY5b3XvV5Hh19O3bl1WrVrld369fP7p27Vrhn9MpCDB//nxWrFjBggULKh1/TEwMoKnww4cPZ/PmzRXOGTt2LEuWLKnV2C83rvgV37NyzYyhCbz43W6s9tqt9gCijvKaiix2dp7IrdJhdyEoKLVzdc8+zJz2OJOmPInDYefXn3/gjrH3VOnY69+/PwsWLGDQoEG6w6wqevXqxWOPPUZWVhaRkZEsWrSIRx55hOTkZJ544glyc3MJCwtjyZIlenmuqpg5cyYvv/wyDz30EO+/rwWIVrfif//998yaNYtffvlFL/zhSVFREaqqEhYWRlFRET/88AMzZswA4ODBg7Rr1w6A5cuX06FDh2rH+Wfgihd8z8o1q1JPYzsPoQdQoU4SdC5VsIzVoZKY1I2bbxnK2Jv706JlS7p3705cs0ZVCv7f/vY37r33Xn17q2fPnlW+T7NmzXjttdcYNGgQUkpuueUWhg0bBsD06dPp1asXMTExdOrUiYgI3zSwOXPmMHHiRKZOncobb7xR7flTpkzBYrEwePBgQHPwffDBB5w6dYpJkyaxcuVKMjIyGD58OKD5NsaOHatXIZo2bRr79+9HURRatWrFBx98AMCZM2dITk7m3LlzKIrCnDlz2LNnj5sz8HLmik/LTUnLZcxHG7HZVYwGwcD2TVhbFnZ7PggBgztG60E9NR1TSWYaEc3iLvpq74rBbiEhLpri4mL69+/P3Llzq3XY1RWFhYWEhoZit9sZPnw4EydO1IXPT+2oSVpuvbDxnfaszSH5cU8Gqiq9hhvW9JY/7MngrhoG9aSk5TLqw98oKQsTvpRMe3IKXbt2pVu3bowcOfKiCT1ogTpdu3YlMTGR1q1b85e//OWivbefeqLq21XpJmTnudi74ZAw49tU2jcNq3blT0nLZeZ3uzlPS8MNtTgf1VqKsUF0ja+d+9nnNIsIqrvB1IA/Y7TblcQVv+I7t+4uJKoqK3Tg8dxJcDoZd6TnV3kvabPgKKg+mlCqDoznTtLamI2wFKCWnKvxuLMLrf403HrKFb/iO0NkZ363u1qhqy2KIogMNusVeSODzcxcsRuLTcWgCGYOS2T3qXws1cQNOApzEMU5OFQQxgCUoFCv56mlRURYM2kZLjEogphwwYlzWZjMQQiDb/3RoWKSkJ/6Q734xru3imTGbQmM+vC3KtXss9/OonH2HzycrHAw38i3sU8Q1Prqau/vUCUvfrcbpMRe5j9wmhN2VfL80l0gqo6zt+dn0tRYQHQThb1ZKqXnMjGZAxEGl69IqqgFZ2luLqRxAwUQZBernDxXlnhrt1Ur+K6+jcspLt/PxeWKV/WddG8VyXUdqraDVUsRccElPN0T/npVKagOn+4tAatdxeqQWvkuDwlXodrwYENIA8wGLR++TaSCgoo9v7zmn7SWouQe56qQIhoHK6hScixP5WieRBoDMTVqgRLgfZ/alUZhAUSGmGkYYtaThGpaesvPn596I/gpabms2Z9Z5TnSYSMiUFsT80slGH1Xm88XYTRzSo2kxCYJMglaRgiktQRHUS6OgmzCSk7SsaFKsElQYpfsPauSVSwxhEZhior1WcXPKrCQU2Qlp8jG8ZOnuGvUGNpf1ZaB1/bmmmuuYdF/apbEUhu2bt2qh/6eL5dTCnBOTg6DBw+mXbt2DB48mNxc71mcU6dOJSEhgY4dO/Loo4/qkYhWq5XJkydz1VVX0aFDBz1K0GKxMGrUKNq2bUuvXr3cOhrVlnoj+N7aZFXAYSciQBP8cxbc1eyLQXAkRwvNOKSkUbBCVJBALcwm1phHfKSCQRFkF2tCX4oZU1QLDKE1iyFw/gVUqTL6zjtI6tGblRu289XKtcx672MOHfW9qUZtSU5O5l//+tcFf5+aUBcpwK+//jrXX389Bw8e5Prrr+f111+vcM5vv/3Ghg0b2LlzJ6mpqWzZsoVffvkF0HICmjRpwoEDB9izZw8DBgwAtEzGyMhIDh06xBNPPMEzzzxzXuOEemLjp6Tlsraa1R60FT+8rH5mvkUiwi78it/l41aVvta67J8rUWX/fGXnJO+CvHnDOkwmE6PGT6TYqpk0MbEtmTj5bwAcO3aMu+++m6IirbLwu+++yzXXXMPatWuZPXs2K1asALTIuOTkZCZMmMC0adNYvnw5RqORIUOGMHv2bP773//y0ksvYTAYiIiIYN26dW732Lx5M48//jglJSUEBQXx2Wef0b59e+bNm8fy5cspLi7m8OHDDB8+vNpIvUudAvztt9/qbb/uuecevfmoK0IISktLsVq1aks2m43oaM0E/fTTT9m3bx+gtSxr1KiRft8XX3wRgDvuuIMpU6ZU2obMV654wXduo/mSiScd9nJV3yJrtOK3bRzCkayiWqf6XkwEgsMH9tGxcxKRwWZKbKX6g9Qg2AxoySo//vgjgYGBHDx4kDFjxlSpUufk5LB06VL27duHEEKvdDNz5kz+7//+j+bNm3utftOhQwfWrVuH0Whk9erVTJ8+XVdxt2/fzh9//EFAQADt27fnkUceoUWLFhXuAZdHCnBGRgbNmjUDtHBlZ3syV/r06cOgQYNo1qwZUkqmTJlCx44d9b/NCy+8wNq1a2nTpg3vvvsu0dHRnDx5Uv/cRqORiIgIsrOz9YmhNlzxgu+M1fcFbcV3qvrSJ7tZoHXZOZ5b8ucQ+rLae2GBRiKCTESFBhBoMvDIlCls2fQbgQEBbNmyBZvNxpQpU9i+fTsGg6HatNnw8HACAwOZNGkSt956K0OHDgW0jLoJEyZw11136am+ruTn53PPPfdw8OBBhBDYbOXVgq+//no9hr9Tp06kpaVVKviXUwpwVRw6dIi9e/eSnp4OaBWC1q1bR6dOnUhPT6dv37689dZbvPXWWzz11FN88cUXXrMZz7eAyxVv4zsDeHz5oNLFxs8vlVCN4P+lawxP3dieO5NbYHfUYTjehURKSm0O2nfoxPY//gC0gpuffvQBa37+mbNnzwLw9ttvEx0dzY4dO9i6dStWq1a/z2g0oqrln7W0tFQ/vnnzZkaOHMmyZcv0JJcPPviAf/zjH5w4cYKuXbuSne0enPTCCy/o2X7fffedfj8oTw2GiunBFT/WpU8Bjo6O5vRprerc6dOnadKkSYX7Ll26lN69exMaGkpoaCg333wzGzduJCoqiuDgYD1f4c4772Tbtm0AxMbGcuKEVsjabreTn5+vT2S15Ypf8bu3imTG0ARmfJuq5aFXhduKX7VzL9RsYM5obY9/4abjKGXVbI1GBYeq4us8sGP8Pmz5GUSaHcRFaNt5pXbJkVyVEmmmoclG60ht+27vWRWLKQxjRM3Dc4VLtd2cYhtXXd2bc4XFzHnnXR5/ZAqA3r0WtJU4NjYWRVGYP38+DofmB2jVqhV79uzBYrFQWlrKTz/9xLXXXkthYSHFxcXccsst9O7dm7Zt2wJw+PBhevXqRa9evfjuu+/0B9j1fZo3bw5oPflqy+WQAnz77bczf/58pk2bxvz58/VMRFdatmzJRx99xLPPPouUkl9++YXHH38cIQS33XYba9eu5brrruOnn36iU6dObvft06cPixcv5rrrrvOv+L7g7DJTHdqKr/2cX42qX2xzkJKWqzfkcOYDtIoM8lHoJY6CLOw5J2kZ6qBNmdc+p0SyJ0tiCWyEqVFL8pQIsopVFCGIb6hAaQFqaaFPn9uJAGIiAokODyQy2KRPAHM+/pK1a3+hdevW9OzZk3vuuUd3Rj300EPMnz+f3r17c+DAAb3XfYsWLbjrrrvo0qUL48aN4+qrtcmvoKCAoUOH0qVLFwYMGMDbb78NaMUqO3fuTGJiIv37969gX0+dOpVnn32Wvn376pNLbXBNAU5KSqJbt24MGzaM5s2b6ynAN9xwQ41TgEtLS5k6dapP50+bNo0ff/yRdu3a8eOPPzJt2jRA276cNGkSoDnn2rRpQ+fOnUlKSiIpKYnbbrsNgFmzZvHiiy/SpUsXvvjiC958800A7rvvPrKzs2nbti1vvfWW192CmuJTWq4Q4ibgn4AB+FhK+brH662AT4HGQA7wVyllelX3vJjdcp0OPltZPfvKFv7jb9/JL2MF/VoZ6f9ZEcdGLKoyKKZfu0a0aBhc45Lb1rPH+HRUe5o3bkB8Q4UQk0CVkhP5kiyrCWODpgij5mRDShw5J+jQwE6QSXC2SCXtnMDUqGWNnI8NQ8zERgb71HzzSqO+pADXaVquEMIAvIfWI68TMEYI0cnjtNlobbS6ADOBymsYXwKc8fqje7akqtoXrgE853zw6v96MIuvNtdM6KWUZK14k2CDnU6NNaEvtUv2ZalkiwhMUS3KhR5ACJSIphzOlZyzSE4XSoSprFFGDXCeHRKglfWODg+sF0IP/hRgb/jyrfcEDkkpjwAIIb5C65e3x+WcTsATZT+vAZbV5SDrgu6tItl4JFtf7aXDjjXjMMJgQhhMSMVQFsATCGiqvlSUavP2a+rJF0IQddMjhJiErtqnFSiI8GYYzd5TZIXRjC0oigM52RjCGmMKrlm9QNdtOqDetdH2pwBXxJdv31tvvF4e5+wARqKZA8OBMCFElJTSzYV7qVto9Y6PQlEEqkOilhaSt+jv/KWDkXMWzYvfqIlCg7IVv9AKWf/vL6AYEAYjwmCi+d/moZgDz3scAc2uIkcN4lieSo4agrFhE1CqVr4MIQ1QAkNqlH3nJDLYVK8E3U/1+PI0+NIb7yngXSHEBGAdcBKosPcipZwLzAXNxq/RSOuA/WcKsJeF7UqHjaahgkUjvdvwZ58Oo9AqyS+V5FskN31ZVLchvOYQ8gKaYfQhscZJbYReEVrKsB8/rvjyJFfbG09KeQoYASCECAVGSikvTPL7efD1lnJ7XDrs2BzwVaoWphsRIGgQKEhoYtDLTYeaBaFmQXPA6hAEKIY6HI3wKZvOV8wGBZvDvZxXaICR6PBA/2rvpwK+PBFbgHZCiNZoK/loYKzrCUKIRkCOlFIFnkXz8F8WuBbHSD1VPhcJReFsSFvGr7MhHXak3UK0zOLkk2GcKZQ0f6uQUDOEBwgiAgVnLQZaXOB2V+eD1WMPUSD8Qu+nUqr16ksp7cAU4P+AvcB/pJS7hRAzhRC3l502ENgvhDgARAOvXKDx1gjnNt6bP+xnxrepbvvrxohomt3zNjET36X5/R/QdOzrevBOvkWihDcm8vHvMDywlIK//oeYh7+4RJ+i5ggEMQ3qRuiXL19e7b7xvHnz9C6yF4LQUO+ViHzBNW33lltuIS8vj2PHjpGYmFjjexUXF3PrrbfSoUMHEhIS9H16TzZv3qxH9SUlJbF06VL9tX/+858kJiaSkJDAnDlz9OO+dvWtK3wK4JFSrpRSXiWlbCOlfKXs2Awp5fKynxdLKduVnTNJSmm5kIP2Fdea+qqsurKutNs8UnJNWgdZgxHFHIghsPYPX10jKHe8ODvcOv+PCjET3ziEqNCAym9QA26//fZKH3AntRH8qsJvLxQrV66kQYMG53WPp556in379vHHH3+wYcMGVq1aVeGcxMREtm7dyvbt2/n+++954IEHsNvtpKam8tFHH7F582Z27NjBihUrOHjwIFCzrr51wRUdueeM0zcIMBqUKre+pWovT8ktrVlm3sXGaFDcbPnIYBPR4YHENwqheWSwTyv9sWPH6NChA5MmTSIxMZFx48axevVq+vbtS7t27fQWUvPmzWPKFC2kd9iwYXz++ecAfPjhh4wbN47FixezdetWxo0bR9euXSkpKXErZLF161YGDhwIaPvpkydPZsiQIYwfPx6Hw8HTTz9Njx496NKlCx9++GGVY5ZS8vTTT5OYmEjnzp35+uuvAVBVlYceeoiEhASGDh3KLbfc4rUrrrcCG0eOHOHqq69my5Yt1f7NgoODGTRoEABms5lu3brpyTae5xmN2ndQWlqqh9fu3buX3r17668PGDBA1waGDBmiX9O7d2+v961LLt+nuw5wBu5sPJLNybwSFlURYSftnsE7vnnQFeG+l9+pWRh7TldM3ayKqnLyz4sXq/avHjp0iP/+97/MnTuXHj16sHDhQtavX8/y5ct59dVXWbbMPRxj7ty59O3bl9atW/Pmm2+yceNGGjZsyLvvvsvs2bNJTq6+dXhKSgrr168nKCiIuXPnEhERwZYtW7BYLPTt25chQ4bQurVnFQKNb775hu3bt7Njxw6ysrLo0aMH/fv3Z8OGDRw7doxdu3aRmZlJx44dmThxYrVj2b9/P6NHj+azzz6ja9eu7N+/v9JuumvXrnXTFvLy8vjuu+947LHHvJ6/adMmJk6cSFpaGl988QVGo5HExESee+45srOzCQoKYuXKlV7/ZlV19a0rrmjBB034u7eKJCUtl8VbT2CtpAqPW2aeRYIPK74C/OMvWsLHqtTT3JzYjPZNw3h+6S72nqmZ8F8KWrdurSesJCQkcP311yOEoHPnzl7LO0VHRzNz5kwGDRrE0qVLa5UhdvvttxMUpAUq/fDDD+zcuVNfnfPz8zl48GClgr9+/XrGjBmDwWAgOjqaAQMGsGXLFtavX8+dd96Joig0bdpUX5Wr4uzZswwbNowlS5aQkJAAQPv27dm+fXu119rtdsaMGcOjjz5KfHy813N69erF7t272bt3L/fccw8333wzHTt25JlnnmHw4MGEhoaSlJSkr/JOquvqW1dc8YLvyp3JLfhvZcLvkplXXYKOK7nFVh4e1JaxvVrq7bpq2pCzsio5AM0bBBFoMujx9Z4jDzAqxPqo3nvimvaqKIr+u6Ioldrgu3btIioqqkqb3jV11zXNFtCTfUBT3d955x1uvPFGn8Zbl6m3ERERtGjRgg0bNuiC7+uKP3nyZNq1a8fjjz9e7ft07NiRkJAQUlNTSU5O5r777uO+++4DtP6BsbGx+rnOrr4//fTTeWHheJIAACAASURBVGffVccVbeM7cXr3F20+XsWK76rq+xYs49kuuyZFP6ojxGykTeNQAk0Giix2mpVl1zXwCMax2FWOZhVdlAq5mzdvZtWqVfzxxx/Mnj2bo0ePAhAWFuZWmSYuLo6UlBSAKttK33jjjbz//vt68Y0DBw7opb680b9/f77++mscDgdnz55l3bp19OzZk2uvvZYlS5agqioZGRl6+auqMJvNLFu2jM8//5yFCxcC5Su+t39OoX/++efJz89388h7cvToUX3iTEtLY//+/cTFxQHoVXmOHz/ON998w5gxY4Dyrr7Lly+vtKtvXVIvBN/Vu18Z0lEL557HpOwZIacAD/aPr3FX3PBAE22aaLsIR7OKyDhXyun8UkICjLRsGEzzBkEEuHQHcu1rf6GwWCzcf//9fPrpp8TExPDmm28yceJEpJRMmDCBBx98UHfu/c///A+PPfYY/fr1w2CoPOhp0qRJdOrUiW7dupGYmKh7vytj+PDhdOnShaSkJK677jreeOMNmjZtysiRI4mNjdXv0atXL59Sb0NCQlixYgVvv/023377bbXnp6en88orr7Bnzx66detG165d+fjjjwFt29PZWnv9+vUkJSXRtWtXhg8fzr///W+9TNbIkSPp1KkTt912G++99x6RkVqx1ClTplBQUMDgwYPp2rUrDz74YLXjOR+u+G65AK+v3MsH645UeU7RvvXMknOYeLWZSctLWFTUg8bDp1d5jQCubdeIx2+4iu6tInlu6S4WuDgQh3SKZu74ZKYv3cWiTceRaJNB59gIHukeTHQLzT70/AaiQsw0jwwm81wpGedKkWXvFeISiVcf02urwpl6m52dTc+ePdmwYQNNmza91MO6qNQkLfeKf1JS0nL5eP3Ras+TDhvhQa7OPRN9lN0cUZuRgXcnlgTWH8xiy7EcFkzq7bVTzntrDpEYE0GAScFmVzEZFWbcloDx3CmiwwMxKIJTeSX6ta6ZdCEBRr1yjgRd2J1C3rpRiN4Cqz4LPcDQoUPJy8vDarXywgsv1DuhrylX/NPi7JZbLR719m5tU8A80yxOyCaMsr5ANt5VRwnY7Cobj2Qzslssi7eewOaQGAyCtfszWb03A7NRYcbQBHKLrfSOj6J7q0j27j1Dk3At0y/QZCC32IoAGgSbdSF2CndGWZcbCeDS784v8OX4Ytf7KeeKt/F7x0dh9MHGdq2wmxxjYEHvwwQIO22VU3xpfpUGuG/PKQKMisAgwGRUdIFeNLkPT93YnruSW2BXtZZaNruqC/3GI9l6B12nmRUSYCQ2Mthr8I1TvXeNzvMLux9PamqyX7FPkDM5p3d8FDOHJfL80l1U5W83hDbkhV0xxBy0MrJNESbFqr/WUTnB5+bX+at1OufQtqO8reLgHjfwzbZ0rDYVIQQFJTbGfaxt9ZmNCovGtCE7O5uoqKhqt278ar2fqpBSkp2dTWCg77UirkjnXkpaLmPm/o7NITEZBIsm9wHgg18O8/vhLAotVRd1NGNjruktBhp2uN9Xbcfd1mcpJpAH+8cz7ZaOldxBY+Gm41p1XylRhFZXT5VgEDD1xnYMaWWqsM/tx09tCAwMJDY2FpPJfRu6Xjn3lmxL1/frrQ7Jkm3pvDq8Mx+NT2bhpuNMX7qryuutmHjA9gSf8QbXGMorjHVXDvKJaTb32p7m9yPZTP58KxnnShnVoyVje1WsKOSs7qtKQEoURSCQmIwKya0b07pVzfre+fFTV1yRgu+pOGcVaMmCKWm55BZbaRJmJrPAWvFCFyyYmWR7is/F6yQr5Z1X+hj2MJe3mHzyKXaka7PrjnRtIvEUfmeSkO7N92Ia+PFzKbhiVf3RZao+gNkgmNi3NR+vP+qbh9+FMIr50vwqSYp7HMCPju78zfYY9rK5s1+7Rnxxn2cpQndfg1/Y/Vxsal1e+89I91aR3JncQl/57Q7J3F+P1FjoAQoIZrx1GntV99V8sCGFOab3MKD5C25ObFbpWB4e1FZ3+L235pDeiMP5sx8/F5srUtUHGNktlm+2pWOza151x3l0tMwnlL9an+Ur8z9op5zUjw81bCIkOITTg972auO74swXsNpVbXtRCOwOzcO/YFJvvzbg56JyRa74UJ6L/+SQ9swclkiAqfoa+VWRTQTjrNM5prr3rRtk+ZmxZ9/W2lJVwZJt6VhsWr6AzSGxleUOOIN/qsNTQ6hOY/BrFH6qwqcV34cWWi2B+UCDsnOmSSlX1vFYa4xzTx2gfdMwNh7JpqDExvzfjlFSiyy6YnMjxlqf4z8BM4kVLpVcUuaBMRBuel3vcONq2wMsTknXw3INipb66nCoevCPJ57Xu8YAzBiawMwVu/XfPTUGV+3Cr1H48Ua1gu/SQmswWqntLUKI5VJK1046z6MV4Xy/rL3WSiDuAoy31jgf/NrkyzsxGASnKRN+80yaCpfVdNMHmvDf8CIpx/MY9/FGLDYVgyK4rkMTvY22AEb1aMmIbrGVOv08BXdkt1g9u9BmV1mVelr/3WpXmbP6gJ4oBO7ZiE6Nwi/4flzxRdXXW2hJKa2As4WWKxIIL/s5Ao+6+5cLS7aln1e+fH6JFi9/XEYzzjqdLOkRv79hDvwyi41HsrHYtBr3dlXy094MjAat9l+ASWFEt1hvt9dxXq9KsJbdx1k70GRUuDmxGWajopf9+vVgFmM+2qir9c5tRM2V4G+o4aciddVC60XgByHEI0AIcIO3G13KFlopabn8d+uJ6k/0kcOyOXfbnmVp8GsE2l1q2619jWHJkreVJH0XQQIDr2pMic2he//HfbyRUpuKImByv3gGJzTVNYDIYLNuFqhAYkwEI8s0BGd/gP7tGnMkq4hDmVrLbKtd5cNfDpPUogG946OYMTShrKS4ZOaK3bRvGlbtqu/feqw/1FULrTHAPCnlm0KIPsAXQojEsgYb5RddwhZaG49k6+2z6oqDtGR29Os8mv53wkWxfjx26+t8lTSN0TuScKhavvzP+zNRVcmWYzmM6BZLqU3706gSPlh3hLm/HkGV2h87NKC8eIUAUk/l67sGY+b+rkcleuYe/bQ3Q88GHNktFrUsndcXdf9y8gv4J6ALjy+qfrUttID7gP8ASCl/BwKBRnUxwLqid3wURsP51TFrHFquMgtAIvjkcAQTrFMpku517JP3vs6CbvswKNpWot0hdZt898n8CqW+9S6+QIFLLoEEvt5ygoWbjrPxSLYelOS8xllj3yC0c512vad54M2B6Io3v8ClwLUJyriPN/p3JS4Qvgi+3kJLCGFGa6G13OOc48D1AEKIjmiCf7YuB3q+eAb11JQGwSaeGNyeQJMmTEqZQEtgm7yK+6xTKZXuCRI9U2dyG7+6HVMl7EzPr9E+qkOVzPg2lchgMyaPyUsCBkVwf794N7s+MSZC3870ZfV27UHgnCguxZbg5TIBXenUVQutvwP3CyF2AIuACfJSxQJXwchusZiMtQtdyCu28cK3qUzoE8foni1xtXaMCqitruEB29+xyHLrSSCZbXqfW5SNbvcqy9mhR1wkZh+1EFWV5BZbWTS5D4M7RbtNYA5VEhZkYkKfOEBzKM5csRtArwGwcNNxr0LsFG7AbaIALsnK620C8lP3XJGx+lXhWv8OnCq77xgEXN8xmh/2ZOjHesZFsu14HnZVcr2SwgemOZhEubpukwb+Znuc1Wp3t3uZDMJNda8Ko0EQFxVCw2ATx7KLKiQZPdg/vkIuQtvGIRzPLdFVf4G2q+DUAKqy699bc4g3f9ivpxE/OaQ9Dw9qW4O/VO3x2/h1R72K1a+KxJgIDIpAERBoUhjWNaZG1zskbD7qrn7mFFl1gftJ7c5jtodxyPI12SQcvGf6J/2V8vx+AT4LfaNQM3aH5FBmIZuP5XrNLFy2/WSFsORDZ4uw2stbZ0s0H4NTfa5Krb6UK69rfoOfC0O9EvyUtFxmrtiNQ9UKY0zoE8f3u8/U+D55Je4loBuGuO+Tr1R783fb31BdhD9A2Jlreoveihb35KuWIYCswqpTiAHOnLP4dE9FCF2IqxJu15Bnf+TflccVm6TjDecKp9nYkt2nz9VJA4ytxyrav8vUawmw25hl+kg/FihsfGL6f9xtfZZt8iqf7l1TQywpNoJd6fluZcac04+iCAZ1aKIfd+0tGBls1ld8zzJifq486pXguxbGMBgUgkwGjIrQttrO476VXfu/xhsItFl5yTRfPxYiLMwzz2Ks9TlSpfe+a+dDRJCpguNCCLihYzRr92fy094MftmfyZ3JLRjRLVYX7MtlD9/PxaFeqfrdW0UyY2gCic0jUKVk9d4MEILRvVry6vDObt1p6oJCi4P5jht5zTbG7Xi4KOEL8+t0EOXNN+rKybruYFaFjkGqhN2n8rE6YwkckoWbjuveev8WWv2jXq34ThvfGUcP4HCoNG8QxPHsIix11PfOkw8dtxGIlSdM5X3kIkUhX5pfZZT1BQ7L5uSt+5ymJ37EEBbFaUMzzIk3ERTXtc7GcDLPvainq6PPs0SYfwvtyqderfiuNr4TZxLLsu0nK72uLvinYwTv229zO9ZInGOB+VVaigzs+RnM7mdhz6iz9LdvxHHuwsc/qRJ2nMgD8Dvy6hn1SvBdvdhGRYt4U2VZsMsFD2cQzLKP5jO7e0vopiKXheZXaGY/RbMw7es4XaBijIj2dpM654c9GYz7WAsw8lYiDCDl6Fl/UY8rjHql6rt6sU/llbBo83FUCRabSoHwrdusQUDfto1YdzCr+pMrIHjJPp4AbIw1/qwfjRVZfH8zhBg1//uZQokhokllN6lzrC5JPJ5BPV93WE/xgV95x/ok7xgD/BrBFUK9WvGhPDgkISZC3+aSQJG16iYbThwStp7Xyid4zj6RJY5r3Y62jYSmoWWCXyQwhjc+j/eoOc6c/XJHn+RRuYikg+/ST+zgHeO/kHar3/F3hVCvVnxAb231360nqG2WbrGPk0RlSBSm2h4gABtDDZv040II7KokomEUQqm8r3xdo0r0nH3NHBJMlV8w0bhKP2ewIYXXxCe0ip9fxZ3+/NSXcOF6JfhONdbVq3+pcGDgcdvDBGBnsCFFP25UBN/fJZlAMQUEX7D3F0CrqGDSsovdcvYfHhDPr51W0Hj/KrfzS5UQOt3yKB2uYGG4nGoSXGjqlarv6dUXaM02QswGHEW52HJPIdXar+a2vDM4Ss75vCdvx8jDtkf5Kc+9Jn/3RlY+M79BMBeur54EbkpoSoBJQUHTNhoGGWD5FBrvX+B2bp4M4a+26RRFa0lGV2oF3/oUz1CvVnzPyL07useSGBPB2z/u59Sun4g/8DmJ0Sa2FkRxpu0IwrreVKP7n136KrbMIwhzEMaIaIwRTbT/w7X/DWW/K4GheodcKybGbe7EmqvT6di4XL1PVg7wiWk2E2xTsXD+NfOCzYYKJkpYkEkv0aVIG+Er/wbK727nZMlw/mqdzkFa6oJwpa6K9SmeoV4JPmg5+bLsfyivfWfLTufOTiamXRvA/6zN4h179YkxnjjyM5hzUwDNw+B04RlOF5ziTKHk9AnJ6T0qpwslp4ol0uQ+MZSeSOUrk42XBrnb9c4+fffb/o4VUyXv6hshZgMlVoeu7ZjLHuyNR7IxSiv/NL7DjYp7mnSGbMBY63Mcpbnb+c5V0WJT+ft/tjO5f5tqG4rUFRfSBnfd9fHb+FcIri2rnTXpnA8xgD0nnY4JmuWz96yKKaF5je6vlhaiWooYHB9Cp8aVO+bsqiSzSHK64CynCzM4U7iT0+GSa1toX8WZQpWmoeUW2ADDTt7OfZb7Tw2DiBgMEdEoASG6xuArZ10y/HrGRdIg2MySbekkRQfwofltBojtbueflI0Ya53OcdmULrERzLgtQRcEs1HR/STHsov17sMXWvgvhg1eXxKT6oXgp6TlMuPbVD1n3mJT+WZbOgkxEahSi5O35aTToZEmcPuyVExRLaq6ZQUcpYWYolpy//9m0iLUQtNQhWZhgmahgmahCk1DBc3CBI2CFWLCBDFhoPUeccdV6J3c2ugUvxnfY8keG6f2S86UmMiUDThraEyWKYbikBiMEdEExHbCEBxR4XpPthzLRQJBlDLU/CbXKLvdXk+TTRhreY6TaFuKe88U6K85V8W//2c7x7LLC4yuSj19wQXf3y+g7qgXgr/xSLZbkQoJfL31BAnNtLLYask5DNZC2jQMQ5WSQ+dMNA6vWa1QU4OmxEz6N+lScrw4H3t+Bvb8TOx5GdjTMrHnZ+DIz0QUZBIdaC2bCBRtYghTmNzdRLNQhX1nHYQHahOG66oe10Dh79e4FvQsKPt3hGKb5IU1Fr40PkdQfHfPoVVAAqEU85n5DXq4tAAHOKTGMM46nQwa6sccDnch694qksn92+grPVTeNLQuqU82+IWmrlpovQ0MKvs1GGgipWxQlwM9H3rHRxFgUvSS1qB10N2Rrgm+LSedtg0VjIrgcI6KPTwGIWq34SGEwBDSAENIAwJi2ld4XUqJWpzPifwMjuZnYjm9n4K1yxjXWbPhb/+qhIM5KjGhgl8nhhAf6T6OXRkO9mapNAsV+uQRahaU2KTPYb4NKGSe+XW6erT+Ph3YhnuKnuasCHcLYTYYKgqZc3VflXqamxObXRQbvz7Z4BeaOmmhJaV8wuX8R4CrL8BYa43zgVmyLZ3FKel6DTontmxXNd+BqWHVnW7OB8+JQTEHUbBlGc3KovZyIjoSO3Ya9vxMRp48xPKgb2gRWKJf3znawM/H7IxabNGPhZq1IJyo/tVH+zXkHItDZhHvOOp2fKfamuVt/sXDLVuyKvU06w9m6XX67uge61XIxvZq6SbwFyP4pb7Y4BcaX1Z8vYUWgBDC2UJrTyXnjwH+p26GV3c4Hxhn++yvtxzHmYVrzzlZLvjZKqaGNXPsnQ/2c5mEmSHELCiySkpDmhIaEokhJJJc2vNX0Z3/yJlEizz9msd6BaBcdR0zjyVq5kR+BpQWoZgCqngnaEwuC82vEu9wz0RMUdsx0T6Ngm35mHfuZsbQBLYcy9FV6pHVtPyCi+N4qy9RdReDumqhBYAQohXQGvi5ktcvWQstJ84J4GyBRa+Ua8tJp2Oii0c/9sKt+J7Y8zPKs/IKK6rrabIpY63P8bX5ZRqJc/rxRyI3YAttxr8cd/r0PjFkscD8Cq2VDLfj2w2dmRYwjfx8TeOw2VVyi601Vqk9HW9LtqXXqZDWp6i6i4EvhqwvLbScjAYWSym9hr9JKedKKZOllMmNG1/cJBRXUtJyWbM/U/9dU/U1D/u+LBXjBVT1PbHnZ+pqfmXpuIdlc/5qnU6uDHU7/qRpMZMN31X7Hi1EBv8JmFlB6Nc5OjO2+EkO5pd/xd7seV9wTXk2KILFKekVavLXNuIvJS2XOasP1JuououBLyu+Ly20nIwGHj7fQV1oXPvoSYcNe94ZOjQKATTBD/ZQ9Wtae78m2PMzaNqwTPALJcZK0nH3yZbcbZ3GQvMrhItym3+6aREWzMx33Oj1unhxioXmV9xbegM/OrrzsO3RCoFBd3TXJr0xH23UVf1F91e/ulaW8uwqpLVZsT3zK5SyieVUXgkpabn+Vb+W1FULLYQQ7YFI4HfP1y43esdH6R117LlnaB4qCQsQnC1SyTM0RDEHuZ1/IRN67PmZuqp/prDqAhypMp4J1mcq9Ol7yTSfUYY1Fc5vL47ztXlmBaFf4ejF32yPVRB6RWgRjR/+cljPabDatZgHX3CmPI/oFluhbHdt4+Bd8ysUoHPzCBCCRZuP+3vrnQd11UILNKfeV5dj6yxXnA6iidfEAXgJ3Ll4jj3VVopanFeu6heCIaxqNXubvIr7bE9X6NP3mvFj/qKs139PEEf5yvwPGrv4BQCWOvrxmG0Kdi/K3vUdtUnnp32Zbsddv9DK1HXX495q8te2QYfrdWaTQmLzCOwOv8p/vvi0jy+lXAms9Dg2w+P3F+tuWBcGVweRE5uHR9/YsNyquZAqPoAjX6ur5xT8DEeoT3n4G9VOTLY9yUemNwkoqxykCMmbpvex2EyckQ2Zb57l1robYJF9EMuaP4UjLb/CPQ0CHhzQho1Hst2yCw2iPK+hMgdbZcddy3j1jo+qscPQOUnPGJpAbrFVnyyWbEv3B/GcJ/Uics+Jq7rpxJadTseWmrDtPeu+ldehaRj7zhS4Cb8i4Npal95yx56vOdsW77Vz4pxkR7Hv0W/r1CQetj3G+y59+gxC8o7pHewYCBQ2t/M/s9/IS/bxGI7n6627nQgBL/+lc4VYfARuArpkW7pua7uGzFYWSutqnxsUwcxhiT7336vKi+86gQD6xOK3932nXuXju6qNTj+2PSedDlGuqn65R39fRkGFFf/2pJgq1UulBrkzga260Gziu/yeOJ3ZygSOt7y9+otcWO2lT59RqBWE/n37bbxkHw8IHBISY8Ldxwy0bxoGlPceUBSBlLD5WC5j5v7Owk3HWZySrv89DIrQW2mfyivBqIgKavzGI9n6RGEva/Xtq01elU/A6UuAS9PR90qgXq34rp7ngxkFLP3jpIeNXx61J9BaWXtyNKuoymaXns0sqkIYzZgbx2FuHFeDT+HOSrU3ATYbb5re9zrpvG0byT8dI3DdlY0IMqGI8rFKiVsjzVN5JaguH8TmkKxKPY3doZlIArgzuQX7zxQw49tUHKrEZFQY1bMFI1268/SOj8KgCD05SpXS58QaX+Ly/Uk7tadeCT6UB/A8t3QXask5QmUhzcPDKbFJjheZaV6WnFOZ/EaHB2IynMNa24J9F4BSzEgvHgmrNPCr2hnPUAxXM0WgOc0ig826am1UBEaXFt4mgyChWTi/H84GtLTmsAAjzy/bpU8edrvK2QILc1YfIKFZOGFBJnrHRzFzWKJbOrRTS6jO1vclLt+ftFN76p3gO8kssGDLSad9mZp/IFtFiSxPzjEoIFXN/kUI1LJV7YEBbXhgQBteX7VXq7Yrz88BGBpgoNBS+3Jfw5T1vGX+AIOXUZiFg3nmNxhjfY7dsnWF1wXQpnEIveKjtBZbztXTIRncKZrGYQGcLdC68H60/ggOzeynS/MIPlp/1E27UYEfyyIhfy2bWAJNmm3+9QN93GxyX/fzq4vL9yft1J56Kfgpabn8uCejQsSea3LO6B4tiWkQ5GavOh+ulLRcdp3M92oK1JQWkcFu+e414S7DGl43foziIvRSlk1WZYSLYr40v8Zo6/Psl+5h0hI4dLaII1lFmo2uCFSHRAI/789k5u2JvLg81U27kWh2vy/oBTxdet2/t+ZQnarn/qSd2lGvnHtOvtmWTsG2FdyWu4CRHbW571COiiE4Eik1OzYhJkJ/YJ3OJOcD5lq5x4kAmoYHENsg0OdxKEBYYO3m3rsNP/CG6SMUUS6UDil40vY35thHuJ0bKQpZYH6VNsJ7mzBVgkOVNAkrDwyyOySzf9h/XiaNwaBwsizCzomv+/lXakHPy4V6KfiZBRYsp/Yzs4+VYR20QJhn+5l5JOR7Sg5rdedyiyuvued8eJ1/PIFmGuQU20jPq1gZNzzQWCHhwSDgH8M7czKvxO14o1AzbRuHVDn++w0reNk0z+2YTRp4xPYIS9V+zLGP5IMq+vR5QyiC0/nuY88pqnndQSexDQJxOFS+8oiwc+4aXNO2ETOGJnhdrZ1beX5v/YWjXgr+rhN52LJP8MVOG+cs2ormUOFUgdS386pyFDlty9G9WmIyaCItEbrX2xWDAqW28iKXTcLM9IiLZHRPTe32FLb4RiFMvDa+km1ByaOGb3jOtNDtqEUaedD2OEebDC47Ini9ij59MVSMQYgOC6gT08VJel4pDqlpE86w3+eW7mLy51t5cXkqGw5lMXPFbq9CXZ/KXF8q6qXgnykoxZZzkhfXWggyattMUW8UsGi3FisfHmj0KSmleYMgHKpmE6tl/7sigMSYCLftv8wCK1uO5bJg03Fe+Da1wvafBF78brd+XAGGdIpGEZKpxq950rTY7fwSaWaS7Sl+VrtzdatIly9U69O30D7I7fxYkcVC8ys0wV3gTueXes/D9ELNynxqxUe+3nqCBZuO88OeDKwOWaVQu2pUzm7GfuqWeif4CzcdR0podNtU2kQqmAyCtDxJgd0ITdojFAPxjapWtZ3oD6jQBNZzxTQZFUb1aImhkqgehyorCJHVrmJz8R9IICO/hBmmL3nI6J4bVSQDmGB9hvVqFwJMCokxEQi39xI8b7+Pbzz69MUpGSwwv0oU5aG7UnqPW/BGj7hIjIpA4B6wJEA/7kqLyCA9G9L1XE8b32nXA3oQkbObsV/dr1vqneB/uv6IVsRStbsl5wS2SKDpX98A4Mw53zrYOFX+vm29F+a8o3ssY3u1ZOawRAyVLJOu4mA0CEb1aInR7WSVuzLeZoLi3tKqQAYx3vos25ROWvnrsnh21UOFUFF42vYAKxzutVPaKSf50vwaERQC2k6A67sqAgZ3iqZj0zC34wYFnrm5IzOHJZYJZvn5ASaFmcMSubad+9/DtRqv83OO6dXSbSvP067ffSofVbprBn6HX91R/7bzyva6bDnpdGxcVnUny+FWfKNlQ9971nVvFcnjN1zF7y45/qA93M7klrG9WtK+aRhLtqWz+WgOhzILKw4LSGgWrv8sAAUHb5jmMtLwq9u5uTKUu63T2EM8ipTsTM9n18ld3J4Ug8moVNhxKO/TZ2OwYZt+vKNynMWh/48Rhc9QIMs/syJgcr94pt3SEdC0pBfKIvQMivY3230q3y3ev3Pz8tr77ZuG6aW7JBWjGSf1ba3f24mnXS/BLTjHNcDIX4Hn/Kl3gj+xb2umL92FLfskHdq4xOi7JOc8c3PHyi73SvdWkYxKbsGCTccBTWhHJbdwezBds9XGfLSxgnBKYGd6PrtPlQkYduaY/s1Qw0a3887KcO62Tmdf2Z68a9jtsu2npdVHAgAAIABJREFU+EvXGFJPnaswudgxMsX2KB/xJv0N5WWx29kP8pn5DcZbp1FMoH7PTzYc5UhWEY3CAsgqsOhCbrWr3PvZZsxGd2WxSXggG49ks/9MAbnFVl0DKSix8cE692q+H68/yuCEpm5/H88ovJHdtPZmziq+ucVWf3huHVLvBH9sr5ZMX7pLS87p5RK80/b8ym2N6Bbrli46opICld1bRfLibQlu4a5OtNVREqTYmWP4J0NcuugCZNKQMdZnOSwrrxnw7fZTlUYSWjAz2fYk88Qb9Fb26seTlQN8bJrNvS59+mwOqdck9ORcqb3CsZ/2ZrB6T0Z5iy6DYNHkPgBsO57rFvTjLWbfuc3nFHTQWndb7SpbjuUwY2iCPzy3Dql3gg/lnXPKVX2VABdVvzarSU3CR3OLrRWE3qhoq3ao0c43Df9Nm3x3oU+XjRhnfY40WXXt/Or8c6UEcJ/1KT43v0535aB+/BrDHj7kbSbbnqxVnz7Pz2N1SD745TC/HjxbQbsxeqnrl5KW6yboI7vFuq3wtSkA6qdy6p1zb+Gm46gl54g2FhEeIMguVsm2mjGUJecIqt7DrwrPCL/K6B0fVWGfvlvLSJ65vgW/Nn+fNvnu6v0xNZq7LDN0oQ/w8BSGmg00DDbRJMy3ba8igrjXOpVdapzb8YGGHbxr+hdGKq7otSElLQeLTa0wKXjzoXja+JkFFhQhUKhZxR4/vlHvBP/TDUexZZ/QPfp7s1SMDcuTc3rE1S7229PjXJUHunuryApbhgZrAZPTniL8jHvJwsMyhlHWGZyi3FNucbhvAxZaHeQU28gs8D3S7hwh3G19lr2qe4/AIYYU3jb9G4WKwUg1JafI5lUDOZRZyOi5vzN96S5S0nIr5PQbFMEvB86iSomiCCb0ieObbemMmfu7P5qvjqiTFlpl59wFvIimbe6QUo6tw3HWGWfLgnecar5nck5NHXtQsVrMjKEJutpamQd64rXxeu+5CAp51/FPOOHevDIrpB2Fty0iaXN+BXs7OjyAM+cs+EJYgIECLxmAeYRxt3U6X5lfpq1SXjj5NsNGrJh4yvYA8gKtDTaHZOGm43y1+bheEcho0HL6BehVepGSj9cf1QOlwO/cqwuq/VZdWmjdDHQCxgghOnmc0w54FugrpUwAHr8AY60TLDa1rHNOxaw8g0epKV/xVFNXpZ6uNuR0bK+WPNg/nq4NrfwY+QZR+e5CT8zVNHr4B+xBjcjzkjfQNDyw0tgAT7wJvZMsIhhnnU6a6l7We6ThV/5h/IwLW3UQPQ1YlWB3qDRvEORWpdcZxOMchbfAHz81x5fpXG+hJaW0As4WWq7cD7wnpcwFkFJmcpliV6Wm6keVq/rOyrpKDXvOO+kdH4XRoGjJOgaFmxObVZuBlpKWy8rft/H/Cp+jSckh9xdb9ILx35JyVjDqw9+8psFuT8+nrmqBZNCQsdbnOKk2dDs+zvgTM4xfcKGF34kqITLY7Fald+awxPIqu0alQuCPn9pRVy20rgIQQmxAMwdelFJ+73mjy6GFliqlh6pfHrxjU2XtmzTI8g319k3DqvVAp+5J5XPlJeI8s+Xi+sGYryAglI1HDmE/f1PbJw6nn+X69UX8PkolMqh8PZho/J5SzLxhH0XNo/RrhgKknsrXi2c66+q1bxrGN9u0en+upb381B5fBN+XFlpGoB0wEK3Tzq9CiEQpZZ7bRVLOBeYCJCcnX5LaVQ67jcDiM8SGh1JqlxzLkzR3Cd6pje248Ug29jIb1KHKCsUnnLiWnLr99DtEegh9fvMBLIp9lfDtOeQWn6GgxL1o5oVAOmzkbViEacdiXrk1gMigilt5DxmXUyLNvOMY4eUOtUegVTI+kFmgVfdRYHFKOnaHVv7rzuQWejzEkm3pWO0qi1PSuaN7bIUJICUtlyXb0skqsNAoLMA/QVRDXbXQSgc2SiltwFEhxH60iWBLnYyyjkhJy8Wee4b2ZRrtgWwVEdpI75xT2608X2q/paRp1WptDonJIPjP+NcJPHeMoLz9AOS1HEy/o+MpPJymt6e+wAss1qzjZK94k4Ghx/jswWBiwhTOWSSzNlj4+zWBNHSpKfJ302JKMfORY2iN38e1nLdBQNMGQZzOK0GV/7+9Mw+Pqjr/+OfcWZKQBAgkIUAgkX0JgoBsKoprQasiVhHUulWhWmvVVqo/qcW1tlZatYh7tYJaRUCKGy64ARIgQAKyEwiQkIQQQraZuff8/rhzJzOTmckkmYQs9/M8PGRm7tw5d+a+57znPe/5vvioD6kaqJru4jjcwb8PNuZxTv8kj1qvw6WxeN0BlmzM89H1v+6lNT6iIe9nHmTx7eNN4w9CpEpoLQUmAQghEtFd/720MBau3oPz2EEGJ3kH9mpG+9NTOzXoRjHmpNeO6V2rpLSxrPfi6j043LJWDlXy3rZyYm5bAV37wdCpLE57lHJXjXJeoN1+kUJKjROZyyj5z2/56xl5fHp9LD3iFb7NdXH6gpM8V3EJN6oPc0L6lhJ7yLaIGy2f1uuzBqfE8+gVGfRLjiMxzo4EDpdUhqVGLNGDsV/+dNTHxZRAlVNj3kc5Hi/KX/nYqUpzY08I6hzxpZQuIYRRQssCvGaU0AIypZTL3a9dLITYBqjA76WULU49YW/hSd/KOX6Vca89s3FxhyVud/QD92gE1KrcY1BUVg1xyXDLp2w4KjmYlY/Vong2qARCSo2Tmz8jdsh5KPbwJb68cZ0opHjlfAZVb+XtW2MYkmTBqUoe/qqav2+JJWHyHLqcNpJs4CbHA7xlf5JYUbNsOM/2b6qx8646KfiHePFTfpln2bIhGNMnA28t4c15pUx/aQ3nDUzGYhE+m6RsFmFu7AlBREpouevl3ev+12IpKXfom3MG1+zKs6XU6OgbRSUaQiDVmEPHKz0uqj9Jbn27DUUWZr6m35yhFhVcJ45S9L/5VB/YguPoPrpePLvebSzf9jUlny/gvlFOHp0Ui90i2F6oMnNJJTs6n0W3m3+NJabmO9goB3Cb835etz3tU6TjSesrVEsbS7WzA32MD5F2WvzP51Qlq7YVYLMILh7SzXOMAN5dfyBg5R+TdpSrvyG3hOIKp9vV90reGaK7+pKGBfYM/Of5CR3s/OOLXQFvfIsikNQE+zxlvfwO7pcUC0KQsymTgv/+iTGJVTz/q1iufm8llXvOJKbv6LDq+6mVZRz7fAFJh7/lg+kxTEzTvYXnfnQw5xsLHc6/l8TB5+o6BX6s0YZyh/NeXgpSp+9jzX+Bp/mR6B1AZu4xFCEoPumo9Z0YlX9MdNpNyq5RDFI9dgiL0N3HncU1JbOEaHiOPuCz9vz2beMoqXAE1OADQEqPCGVCBzt2qxIwjje2T1dOS4zFlpSGpUNH7h1vZ3QPC29OjaHkk/molScCvMuXyn2bOPLanVxlW8PmWXFMTLNypEzjZ/8p5w/bBpBw4/PEDjkvoNEbrNaGc5fzbpyypqCnRUj+aXue85WNQd/XHBhJTBI9RbgogNEDhHSn2iHtxvDLKp1oFaWoVScZ9EI5sU+UUSmjsMTrOfB3nNOn0W6g9yYdbxlpf+ktQ4TS6dLIPlzKsJ6diIu21DL+FVsO88X2AhR7DImX3sevV1ZzpExjYpqV3w2v4Ngnz6MFiQBqzmqOrVqIc/nDvH1JJW9NjaFTtGDJdifDXqrmx7RbSL72Mawdk8K6ts+10dzjV6fPJlQW2OZzttLwOXxjCTeJSVVN0U5v2o3hr9lbjFpRirVTN0BQreKzOcdfEaYxeJd3nj6md63wvKdop4DF6w6wfn8JZVVqrZGqtNLliX5H9RyEa8Q13LxMl+N+dFIUA06spTznq1qfX52/myNv/JYxx1ayZXYc12bYKKuW3Lyskulfp2C/Zj4dz7zCc+3h8j9tHL933oHmZfxRwsXLtmcYK7aHeGf9OP7t2xStnE/5tq9RyxsejR+cEo9VMeTPTVffm3Yzx3e4NOxJafSc9SqasxrX8SNIh66tZ61Pids68F5TVoA+yXG1RyUh9GlHPbPyOk2Yzlf/yeS5H/P4zRg7b18Vw+g3XiS6VwbWTslITeXE2vepXLuIv0yycu94fQfgDwdd3PBhFUUDr6L7jTMRlvrvtzdYok0kyuXkSdurnudihINX7X/lRsccNsoBDT63Qd/CLxAni8jetooiFbrNeIroXhn1Pk9qlw7sPHoSiUSVsCO/zAzuuWk3I77Ta0lIsUVhT0onqucgADrHNtwQvNmQW8K8j3I8iSQaBNTX8xfEDBdhsZJ42X3M+Upje6HKkCQLT01UKfrf33EcO0TB2w/Qe+fbrL81mnvHR+HSJA9/VcWkJXGUT3mShHNvapTRGyxWL+DPzht8nosTVbxhf5qhYl+jzu0qLeD+YWVk3h7HjGE2UKzYu/Vr0Lk+31bgWQpUNcnD9SjTDW27mk+7MfyTIdJfT+sanpx2KIytuVvySkMeZxH4qejWD1vXXkSdfTMzl1TiVCV3j43iXPt28l+9k7vS9vLjbbEM62ZhZ7HKhFfLmX/8PLrd9DzRqUMb/JkG3o7R6+pk/uKc7vO6UadvgDhIQ6ncn8WIFP22zMpXieo5qME5C/6ommTJxrxazwcy8LZezaddGP6idQfILwu+d31At4av3xsYy3J1jeXGNtTGED/yUrZHD2fu1/o1vTk1htU32vnbxdFEWQULMh2MesvC/nEP0nXy3ShR4asGh8LfUVmgXs4/Atbpe5zUojWUZX2C83h+vT5DObiBgYkKTlWSU6gRnT6isc324esdRz0CIBDcwNt6NZ92Mcf/OPtIyNeH9ugU9rmC1Xb3Xse3KIIRvToH3E4biYQWIRQ6DD6Hpz/J4pfDbQxKtNAtVqHgpMaty6v4UjuDrjfcjSW26eezz7qmEYWDWdYVnueSxAne7bqQcz4oZt+nEmvn7kSnDyc6bQTRaadjiekY8FxSavSvzEYRgq2FKg4VuqRF1vAPHa9i0boDvJ95kEcuz+Dj7COeJJ9qp5516b0q01bFPduF4U/O6O6p2R6IUAUyvfFX2vFOAfUX2wS4ZuEan3TTSKBVV3Dsi5ex7vycN6+M8QiKADz2TTVfuE4n+RdzQ67LRxbBU67riMbBTdbPPM/2jHHw5S9jmfh6OQePH+Fk1hFOZn0CCOwpffVOIH0E0alDEFa3su/R/ZyeUAnEkJWvIqJisXfv3yStdqqSucuyfUQ+JPruwIwenXwkwtuiuGe7MPwZY3vz0Idbg4624fbmgdw/f4lo47jNB49H3Oir8nIoXvF3zupcxJuz4+jdSaHcIXk3x8ktZ9h54oJo/vdiFhV5OQ2KgjccvU5fFE6us9YsL6Z3VthxVxyLtjr5/qBKVr5KTqGGI383jvzdnFj3PsJqJ6rnEKLTR6CePMaIFL0jyyrQiO49AqFYgn1oo5Dooix662s8MZdL83QIbTm/v10Y/qJ1B4IafZcOtrB/WB933qv2u38ZKIdLi+iOWqk6Of7d21St/4DHJ9m5f0IHFCFYl6dy/dIq9pRFEWd3cs1QG29OjeaCJc/Q7eYXIja3D6uNKDzkupVo4WCq5XvP8zE2wa0j7dw6Un/sVCXbizSy8lU25ev/b87PoiQ3C4ARt+ptzspXsffs0yxttyigKAqqqoEQng6hLef3C9lUez/rYPTo0TIzM7PJP2dDbgnXLlzj+TH9UQT8d9aEsH/cDbklLNmYx38zD+LSpEdcs6TCwaHjlbxjiERGCEdhLkUrnmGA3M9/rophRIoFlyZ57BsHf8lJpPOlv0erLse59CG2zo6jZ0eFh76s4h8lE0m8tPmlDy2o/NP2HJdafvR5/rjTytFyjX6d1IASZ7nH9U5gSn8rNotgxItlbC6QWDt106cERnygQ/jxmHC5eEg37ji3L0s25vHu+gMe1SO7VWHxr1r3iC+E2CClHO3/fJsf8Q11nGBIWb/NOaPSEjzn1KQeEHp4qT6NsCgiYnvopdQoy1zO8dX/5u7RgqcujCXaKth9TOP6JZVsS5lC0o03odj0XX6VGVdw07IVfH5DLI+cG8Wnr37JgR1j6DBwQmQaFCZ6nb67sDPfp05fZ5uL/R37cGbJbSQe/ZGBlVkMUQ4wvIuD07tZSOuskNa5ZpEpa1Y8x6skm/PLyCr4mqwDX5L1o8pukY7SawTR6WcQlTrEc/2N4YufCkiMj0JQs3Ih0IuetmajD0W7GPGnLfgh6OsN6dUDKb5EEn3P/LN0ObaVN66I4aK+ev/88kYH93/XgeiL7yXGb5lLuhwceeMenh5dwD3jothRpDL6LQudb3gBS1z41ybcE97GXlkUDl62+dbpA1ivDeCX7jp9UkpcJYdx5G6i256l/KpfCb8dF0XBSX3I7RZXe7XZoUq2FereQdZR2Kb1YmfsKKpTx2Dv1rdRMQG7RYAQqKo+lQsk8dXaCDbit2nDN5be5q/aGXTt/Impw5gxtv4CHA99uNVTJDNSSCmp2L6a4s8WMK1vNQsvi6FLjKCwXOO2j6pYpUygy8W/xhIdF/D9joI9lCy6j/W3RpGRbOFf6x3c/1MGSVf/KawovwD6JsXSJymOLXnHw9btD0Y01fzb/hfGKj/5PP+DOsSnTh/AkTfv5c9DcplzdhR/Xl3NY1uT6eo8wogkyYgUi/ufwoCuSsCpwspdTi5fHkPP2a83yvgvGtKN5Pgoz1TOW/uvNXYAwQy/zSbweCdmhEqYCXcpz5+rRqZGNMdfrSyjaPnTOD77G6//TOW/v+hAlxjByl1Ohr8m+WHA70i6/A9BjR7A3q0v0eNnMnNJJQ5V8usz7YyXmziZ9XFYbZDA7sJyPttWQLmj8WW0qojiFsfv2aj5ptxOsGxjoe1Z7OjZlGrVSRz5u2si+vkqKdOfwH7bO2wa9ScWOC/jl9/1YPAL5cQ/Wca4V8qZ/b9KFm5wsC5PpdIpKTgpsXXt1ehVgNU7Cz0Rf03WaP+1tey9NjvH9xG4CEJjdmyNSktg3hUZzF2WHTKGEA6V+zZRvPJZxncp5a1ZcaR3VqhwSu7/rIrXCwfT9drfERvm9tmOY65ix+71/PGL3aTECb47pBDXr/4a3cJvXUIRtTP3wqGcGG5yPMDb9scZpuz3PH+eZTPP8Rx3Ou+m4sAWkJonVTenOsUzPYk5bSQxp+lLAmpFKVW5W8jJzSJzVxZqpq5SbBEQa4focY1P9nGpGjmHSrEqAqdbIzGQgk+wRK7WQkRKaAkhbgL+ChxyP/W8lPKVCLYzJIF+BGPprcoZ/Ka/YFByg3+0DbkllFQ4mHdFBjmHS1m87kC9q81pziqOr/43VZs+4tFJUTxwlr5Ml3lY5fplTo4Ou4nk839er+2zQrHQ9bJ7efb132DtnELi9fdhT0oPeGwoYz5R5UIRkBwfxZUjenLR0BTW7i2mrNLJS9/urVcnYNTpe8f+GIOUmjz+SyyZ/J0F3LDfQkqcICVOobRKkt91JF0CnMfSoROxg88hdvA5enzgeD5VuVlU7c/iZO5mkhuZ5WcE97bklWKzKlw4JJnVOwtRVd/svVCJXK2FOg3fq4TWRegy2uuFEMullNv8Dn1XSnlXE7QxJMF+BCOT7tZ/r+d4ReANOnec2zfoOdfuLSahgz1g5lagWnmKRaDVI9hXfWQXRSueoZ9ymLdvi2VkdwuqJnns22qe+qknHX9+Px0TGyb+aeucQrfpj2NPOg1hDbwbz6IIfn56d5ZtPhx0JUKTcMHgbsyZMphF6w7wWU4+Ww+VNmjkP0481zse5F37PPoqNSnUl1vW8Fx/K++U1bj50eln1Hk+IQS2hO7YEroTP2IyUlMbpbJjswiGdO/IlrxSXeBT1RjRqzOzzu1ba1CpK5GrNRDOiO8poQUghDBKaPkb/ikh1I8wKi2B9C4dyKqovWPuyhE9Av5Y3katSX0UiLLpxp1zWL8pdhWUeTwJp0vjte/2+ii8hkJqKqVr/0vp94v59SiFv10US4xNsLdE44alVWT3vIquM69r9PbZqO7B98X3S47jlrNOY+6y7DqXH7MPlXLPO5tYmuVfSqH+FNGJGY6HeM8+jzSlpsrajYNdDHCX+M4qkEQPrn/WYWPn9n++PIOBKfGe2gdQU87L/z5pC3n8kSqhBTBNCDER2An8TkpZa29mU5TQqutHyD9RVes9g1PimT898KjiHxswijg8vCw7YAquBPYUlofVVvVkCUc/fIwuJ3ayeHoMk/vrX//rmxzct7YT9ovmkpA6pI6zNJ4LByVTUuEIKNtls+giIUYSy5a8UjbXsdW4PhTQhZnOh3jXPo+eombH27hU/bvIqU72FDhpTnIOlzIwJR7VvZSpSnhkeTYDU+JrGb7/vozWNtpDeFH9cEpofQSkSylPB1YB/w50IinlS1LK0VLK0UlJ4QWr6sJf5NL/RwjU+FAy2uP6dK2lkSeECJp3r9VjzVuJjmNgXAVbZ8cyub+V4gqNae9VcPf+icTNeI7oZjB6gFe+2xdQ5FMAvxjdi3fvmMA5/RNRRNOUy8yTScx0PEiB7FzrtTP7dKS5inR6I9HrInj/zkZRjkB46yu2RsIx/DpLaEkpi6WUxqLvy8CoyDQvPEL9CNUBilksyzrMohBr8N4KOQI4f1Byrc6gIQirjZKz72dHseTT3S6Gv6Hw7eA/RnTPfDhoUnp2n52e2gmbRWAR+pTGSFi558IB2K2K5wYR6FNohchU9tovuzPT8SDF0rcTLk8YFIGz+5LWpQOhfj6LgIwenWp1N0ojlZdbMuG4+p4SWuhR++nADO8DhBDdpZRGxOZyIHLKi41gQ24JpVW116MlMHdZYDdu7d5in+CVRRHMOrcvkwYm839LtzY6D9+S1JfrN07l6L6ddLn2bjrE1h71mhLFXW46oYOdeStycLj0ApXTx/T2SVLxdme9g5wA81ft5PvdRY3+LnbLVK53PMhi+2N0FuU87byWf1n8K7A3ntxjFYDeYY1OT2BzXimqqqEoAs29Xj9vRQ5zLxuK3aIv4ykCHr1yWKsd0esiUiW07hZCXA64gGPATU3Y5rBZu7c4uIvurmobKHATZVNwOPUbY94VGZ5jjLVdAbWW7hQB3TtFU3TSEdDL8PnsM64h8QzRjHvmdQQwfUxvMnp04uPsI55YhqpJenSOCTiXDXTj33PhANbtO+ZTGsyi6Pse6tsZbJdp3OiYw3BlD2+pFzfkssJGAv27xTNn8mBPpaN3fjzgEeHIPlzK4tvHt+q5e7hEqoTWH4E/RrZpjWdcn65BK81Yg0RjAwlqvPDVbg4fr/SUwhYCrO55v3FuTerqLuFQX1nrurAokNIxOqzPLyqr9kk6EuCp/GPUpa/rhh+VlsDVo1JZ7N7uLNBd5a2HGhYE3CL7skUNvLQaaQQ1uglLNuZ5chkk8F7mQaaNTOXOSQ0T92xNtNnMPdB/4EEp8T6lmD2EWMcyRrpF6w4w1x3Nt1kEVovi2cBx3oAkJHqKZ6CimM2JpsGQHp18DH9MegJZeaU+bVMUwaptBT7eihBw0/h0j9sfTkLKhtwSvcOw6J2fzaowvk/XOqP/Nouga6y90XsAGooidJk179/V2+lyqdIjvQW+iWFAm/IE2rThA/Tq0iGg4btUyfxVO7nnwgFB1/O9R0aHKrloSDIC+GJ7AZ9vKyDKpvDIz4fyzy92nrKbGfTRqrTC4fFuBHDuwGQecLu0ZZVOXvluX9DU4pwjJ8JOSPHOc7BaFK4d04tpI1NZu7e4zjp+LlXSrWN0wO/KbhEBdzuGUxsw3PdoEuYuz0Zq0lPrwL//N/oBn+tU9F17LrX1Zur502Y36dSFBny/uyjo5gs9yOd7V3z5UwGfby/wrPU63PPC3l2aLyIfjGPlNTXjJHrJMGO1Iz7GFnDNXqAH+iZndPeU+6orIcU7z0FVNXq6YwNGbMQiINqmMGtiH4andsJbSdxiEVx7Zm99+6sfwbY4NyQMUlfn4xO8FbonYnwXV41MrXWdTlXibGOKu21+xE+M9xVqGJwST2J8lCcqHWyEMxKDvMtc+1e+kcB/Mw82Wi67LhShj0zS7znvG7hPUhx7Css9x7z83T56d42lpMLhWbM3JMOuHpXqEZQ0XNeBKfFhubLBEqaCJbUsWnfAk/ykCMHAlHgW3z6eX725nmPlwWsdGERYthCrRaCg775ThB68DXTt/qrJxj791pqp50+bN/yMHvqoo0rdnXxs6jAA1u8/FjLl0riRDZmtQCOSpPEa+XXRMdrK6zeP8bRDn5cKnxHckI764qejtSrHSOkrDxbMsINF8AMdFyxrLdA5SiocGJoPTpcuX/3E1GH0S4rjx0bUxQsXIz3ZSLee5jWie7e/ruy8QO9pzbRpwzeCOJrUl+IeubxmaS6clMtRaQnsyC8ju3spReUODpVURqxtdqvCLRPSefGbvSGPO9+9g3BUWgJXuefSxhq80XHdcW5fRqUlcP6gZD7fVuB5r9EJOF0aJRWOiEWrw+0kQB85rYo+f5fUyFdnHTxe788VQEKsLainMCY9gWqX5hNkHHNal4BCK+F2coFUlNsCbdbw/YNzRraaQTg371Mrt/sYZkMCTUGRkrLq2slFdouCw2tOsTI7nxvcSr7ebTbc04QOds+cc9a5fVm94yhOVfoox55K93RUWgK/GN3Lo3SsqhofZx/xWU680K16U1hWzdfu9vt/zwK4bqyeg/Dgh4HLcle7tFqrCxn1KJbSnmizhu+fvKOI+olubMgt4aVvfUfjSDr1qiY5GqCsl8MvkBAsBmE89t8e/IvRvUK6tKeCq0am8sHGPI+HMjmju89Ua5bbY9mQW0KSVwfgUiUaNdmGxjXZvLLrvGMfWw+VknP4hKeDVkTDFZbaOm3G8P3FOIwoc7V7++zI3rVTY0OpqPin7tZFXd6AQNdzM25ogK92HA3xDh0lhEqQd+S52qnxf0u3IqVvzn2gQrItAAAO5klEQVRLcE8DxQX8A2r+GgePXJ7hCUx6xyZe+Gp3TYcu9e9H9Xh1+pO62rGs09Np7So6jaFNGH4wMY65lw31RJR/3K8r4y6+fbznRrvu5bWeUcdfaXdcn64hFWq8X7NbBFOGdQ+5Z11RICk+ilvOOq1mTb2OwKBFETx6hb43/YWvdtcyAiPy7HBqaNSsSTucLU8cItB82fuxv65CsJiEd7TdP8hpLMkFCmT6G3lbUNFpDG3C8IOJceQcLg24zXJUWgJLNuZ5stocLo2Fq/cwvFdnz40xKi2B28/pEzD4Zhikd6S4rrVdVdODjRal9hZfhdq5/xYBj7qXmgIJgxg36tu3jWP+qp18t6vI43GE8hL8aSmjXrjiFv6bh4wgp78ctnFdBv5G3hZUdBpDmzD8QDfNonUHeOdH3623NkuNQfiPtZ9vL2DV9gKf3n/OlMEALPxmr8/xUtMDhY+7lwYBduSXBTRqbyQEzJ67cEg3vtlVWEsfsKTCEVgYxKn5ZB3ec+EA1u/XN80Ya9Ph3MQtadSrj7hFoCCn/+jufV3TRqbWMvJxfbpitSj68+gJT+2JNmH4gdZc5y7L9vGkh6d2Yu7Ph3pujmkjU3nfnXwjqQkSOfx6/zlTBnPR0BQWrt7DF9sLkBLsNt8RaUNuCfNW5KBp+vwypWNUWBtmjNH7jnP7eko4GWv13qOe3ar4jPhG1uH6/cd8Rv76jtynatQL5mU0JCYR6D3+1yUhoDehanpHKyUez87o7Ns6bcLwwfcGePDDrT4jq1URPkZvHL/49vHMX7XTp4R2oOj/qLQE7ji3L0nxUR7X3vtcH2zM82T4SSnJLw3P6M/un+izV8B7rd7bKLxd24+zjwTMOmyI0ZwK7bjm8DL8r2vayFTPdMw3SOj7vhe/2cvi9QdIiovilrP7NKjQSmuhzRi+wYbcEt7fkOd5bBEEdX193GSv/feAzxbVQK6j8VkfuEdpo5vxrr8WCotFBNwgFMiAjcdr9xbXWgprjLE21FNoDM3hZQS7Lv/grdUiaomklla6KK108eCHW1m4eg9/v3ZEm5z7tznDX7u3GJe7KzeEJ0L13Eb0/+PsI0zO6K4rrfpF+/2XzT7YqHcsM19Z65PLL9DlqL/ZVejJ8Q6U1iuA8wcme4JP3jdWIDc4kJx3qPTb+tDcS37N5WXUdV2j0hJ49/bxIWsg5h6rYNqCH7hyRI+g4qytlVZt+KEKaRg31lUjUz2lrQO56cb8vNqpsWZPMSN7d/aJ9i/ZmOcpl+Wddmq87m303vN1o12f5+SzNOsQCR3s7Ckq9+zn/3rH0VrBRMPAq516pzHvigxmjO0d9lJXqO+lpXAqvIxQbbnlrNPqTJtemnWYpVmH21QH0GoNv65CGt6BvuteXusx5vczD3rW8kH3EIxR26VJ1vtt0ZXoN8h5A5P5zJ0Hr6oauwvKfFz6i9wbZbzdykXrDnhuqvwT1cya2If4GJtH8snf3fVvi6ELWJ9RsiVF6oPRUhKLQA/mrd1bTFYYEuJLsw7zWU4+b7XA77S+hLUfXwjxMyHEDiHEbiHEnBDHXS2EkEKIWtU5I02guaLBqLQEz4g476McHxUaf8nkWnLaUo8LCPTEnGluj+Fr7yw7IVi/v6aDEMDwXp1r3QwfZx/xeZxz5AR3TurHtJGpAfe/+7dFkzV5B6EkxMP9XkwCs/Sus7lyRI86jyv55k2Kd6xn2oIfQqo0twYiVkJLCBEP3A2sa4qG+lPXKBishr33Wj74Fr/UgmxhfeGr3T6rBJrmu4kkkAzzhtwSYmy+1V0mZ3T3fGaw4JN/W7z3u4czyrSFKi+ngvnTz+CG8eksXL3H49l5U7lvEyfWvMcJIDbjAh6ovI3HVuS02tE/kiW0HgWeBu6PaAuDUNdcce3e4lpBNf+1fIMZY3uHFKJI6GD3cev9Q0GKn2i7t7tt1GS79kzfIGMwQ66rLXXRkubQrY1RaQm8dONoFq074LMDUHNUUvzJc1w2wMo1Q6z8YdWXHH51A10ums20BRojUjux9K6zT2HL609ESmgJIc4AekkpVwghghp+pEtohRoFx/Xpis1Lx81uVQIafTjnyj5cM/8T4MnQMzoAf6lub3dbaJKLh6bUa024sXPgljSHbo0Yv9VDH25FAlW5m1HLinhqZgxDky1cMcjGI19X8NzyJ6noN4ENF80mfU4pT0wd1mrW/htdQkvoWtHPAvfVdaKmKKEVDCNBZ8bY3swc27vWJpxw8c8LsFn0aPt1Y3sH1akz3G1dsUmQ0MEekWsyaT5mjO3NvqcuZWL/RDr0H0fK9X/lqk8TWL7DSccowd8viSbrjljGONZx+NXZlGV9yj3PvMHEue8E1HBsaQhZR7lUIcR44BEp5SXux38EkFI+6X7cCdgDnHS/JQW9qMblUsrMYOcdPXq0zMwM+nKL4YWvdvPMZzs86bIzxvb25OiHWjbzlnD23lRj0vowXH+pOild8x5nFX7APy6x06+LPm6+m+3kvs+qsCqw7a6OrFDHoV6xgOvGnXaKWw5CiA1SylrB9kaX0JJSlgKJXh/0NXB/KKNvTQTKCzAI5VIb1Wgl7XP3V1tixlg9PrNo3QEetNjYVHgWIz+cz29OO8BD50RxbYaNywZY2X9co4NV4trxHfe+8AF/+mgIV4/uVSt3pCVQp6svpXQBRgmt7cB7Rgktd9msNk19ltK8MTqMcCSrTVoHM8b2Zv9Tl2JPSifhumf4l+V6Bi2o5v1tTmLtgqHJFjRN8tFP1ZR8+QrlhXksWneAaxb+0OLc/zpd/aaitbj6jaElZ9CZNI4BD63EoUpObFhO+VcvsefuOFI71oyjH2538rtVKqXDZtBxzFTiou3kzPtZs7czmKvfbgtqNAdGIpFp9G2PnY9P0TMxR/6cP958MakdFX4qUpmzqoqTDsnUwTa2z4riN3IRJYvupSRvF+lz/tdiEn/MEd/EpJFc/Oj73FLxKm+vK2TVXhddjm3hrxdFM2OYDYA9xzTu+czBt12m0mnCdHolxvPdnAuapW3BRnzT8E1MIoAn8i8l5du+pmTVS5zdrYLnJ0czrJuewblip5N71nSl4vw52JPS+WD2hCb3Bk1X38SkCTECfwOS44gbOokety0gM2YCI18q57efVFFaJZmUbsV58hhKdBwA0xb8wI2vNkuGey3MEd/EJMJsyC3h6gU/IIGKnWs49tm/6MpxzuxhYU2fWcSP8A3y9U+K5fP7zmuStpgjvolJMzEqLYF9T11K/6RYOgwYT/fbFlDe5wK+cAwlbvgltY7fVVjOwIdWNuuSnznim5g0IU+t3M6r3+/TRV1dToTVFvL4SIt9mCO+ickpYM6UwexyL/3VZfSgi308tXJ7k7fLNHwTk2ZgzpTBetafJdCeN18WfrOXK57/rknX/E3DNzFpRnY+PoUnpg4jyhrc9CSwOa+UBz/cyuT53zTJ3N80fBOTZmbG2N7seGwysyb2qfPY7fllXPPiDzz04daIdgCm4ZuYnCLmTBnMB7MnMGNs74CiFwaq1BOEZr6yNmLGbxq+ickpZFRaAk9MHcb7sycwJj14Fp93zcRIGL9p+CYmLYBRaQm8N2sCT0wdRr+kWJLj7fTsHM2Z6Qm6mhM1NRMjMfK36HV8c1uriYm76MtHOWx2a/8rAs7qlxiwBJs/jVHgOSW0hsIQJibNxfYjJzx/axK+2+VbLbm+tFhX3ywMYWKis3ZvsU9dB8BH0q0htFjDN6WrTEx0vG3BblWwW0Sj7SKsOb4Q4mfAPwAL8IqU8im/12cBdwIqutru7f6Vdvwx5/gmJuHjbQtA2HbRYCEOdwmtnXiV0AKu8zZsIURHKeUJ99+XA7+WUoYUGDM36ZiYND2N2aTjKaElpXQARgktD4bRu4mldpUpExOTFkRESmgBCCHuBO4F7MD5gU4U6RJaJiYmDaPRJbQ8T0j5gpSyL/AA8H+BTtScJbRMTEyCE47h5wG9vB6nAodDHP8OcGVjGmViYtK0hGP4nhJaQgg7egmt5d4HCCH6ez28FNgVuSaamJhEmjrn+FJKlxDCKKFlAV4zSmgBmVLK5cBdQogLASdQAvyyKRttYmLSOE5Zrr4QohDIbeaPTQSKmvkzmwPzuloXzXldaVLKWgG1U2b4pwIhRGagNc3WjnldrYuWcF0tNmXXxMSk6TAN38SkHdLeDP+lU92AJsK8rtbFKb+udjXHNzEx0WlvI76JiQmm4ZuYtEvapOELIX4mhNghhNgthJgT4rirhRBSCNEqlozCuS4hxDVCiG1CiBwhxKLmbmNDqOu6hBC9hRBfCSE2CSG2CCGmnIp21hchxGtCiKNCiOwgrwshxD/d171FCDGy2RonpWxT/9CzC/cAfdB3Cm4GhgQ4Lh74BlgLjD7V7Y7EdQH9gU1Agvtx8qlud4Su6yVgtvvvIcD+U93uMK9tIjASyA7y+hTgY/SNcOOAdc3VtrY44tepH+DmUeBpoKo5G9cIwrmuXwEvSClLAKSUR5u5jQ0hnOuSQEf3350IvUmsxSCl/AY4FuKQK4A3pc5aoLMQontztK0tGn4g/YCe3gcIIc4AekkpVzRnwxpJndcFDAAGCCG+F0KsdUumtXTCua5HgOuFEHnASuA3zdO0Jieca28S2qLhh9QPEEIowLPAfc3WosgQji6CFd3dPw+4DnhFCNG5idvVWMK5ruuAN6SUqeju8Vvu37G1E5bWRVPQFr48f+rSD4gHMoCvhRD70edWy1tBgC8cXYQ8YJmU0iml3AfsQO8IWjLhXNetwHsAUso1QDT6RpfWTn21LiJGWzT8kPoBUspSKWWilDJdSpmOHty7XErZ0pU/69RFAJYCkwCEEInorv/eZm1l/Qnnug4AFwAIIQajG35hs7ayaVgO3OiO7o8DSqWUR5rjg1tsJZ2GIsPTD2h1hHldnwIXCyG2oUud/15K2aIrkYR5XfcBLwshfofuCt8k3WHxlowQYjH6tCvRHZ/4E2ADkFK+iB6vmALsBiqAm5utba3g+zMxMYkwbdHVNzExqQPT8E1M2iGm4ZuYtENMwzcxaYeYhm9i0g4xDd/EpB1iGr6JSTvk/wE701Ni8BKEZAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "(handle0,) = pl.plot(x[:, 0], x[:, 1], \".\", zorder=-1, label=\"empirical\")\n",
    "\n",
    "for c in range(num_clusters):\n",
    "    (handle1,) = plot_cov(\n",
    "        gnd_mu_[c],\n",
    "        gnd_cov_L_[c] @ gnd_cov_L_[c].T,\n",
    "        color=\"k\",\n",
    "        lw=4,\n",
    "        label=\"gnd\",\n",
    "    )\n",
    "\n",
    "(handle2,) = plot_cov(g_mean, g_cov, color=\"C1\", lw=4, label=\"Gaussian\")\n",
    "\n",
    "for c in range(num_clusters):\n",
    "    (handle3,) = plot_cov(\n",
    "        mu_[c], cov_[c], color=\"C1\", lw=2, label=\"EM full rank\"\n",
    "    )\n",
    "\n",
    "pl.legend(\n",
    "    [handle0, handle1, handle2, handle3],\n",
    "    [\n",
    "        \"empirical\",\n",
    "        \"gnd loglik={:.3f}\".format(\n",
    "            mixture_full_log_pdf(x, gnd_mu_, np.linalg.inv(gnd_cov_L_)).mean()\n",
    "        ),\n",
    "        \"Gaussian loglik={:.3f}\".format(\n",
    "            Gaussian_log_pdf(\n",
    "                x, g_mean, np.linalg.inv(np.linalg.cholesky(g_cov))\n",
    "            ).mean()\n",
    "        ),\n",
    "        \"mixture loglik={:.3f}\".format(\n",
    "            mixture_full_log_pdf(x, mu_, kR_).mean()\n",
    "        ),\n",
    "    ],\n",
    "    loc=\"upper right\",\n",
    ")\n",
    "pl.axis(\"square\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Environment (conda_mxnet_p36)",
   "language": "python",
   "name": "conda_mxnet_p36"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}
